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Abstract. The spallation reactions are a type of nuclear reaction which occur in space by interaction 
of the cosmic rays with interstellar bodies. The first spallation reactions induced with an accelerator 
took place in 1947 at the Berkeley cyclotron (University of California) with 200 MeV deuterons and 400 
MeV alpha beams. They highlighted the multiple emission of neutrons and charged particles and the 
production of a large number of residual nuclei far different from the target nuclei. The same year R. 
Serber describes the reaction in two steps: a first and fast one with high-energy particle emission leading 
to an excited remnant nucleus, and a second one, much slower, the de-excitation of the remnant. In 2010 
IAEA organized a worskhop to present the results of the most widely used spallation codes within a 
benchmark of spallation models. If one of the goals was to understand the deficiencies, if any, in each 
code, one remarkable outcome points out the overall high-quality level of some models and so the great 
improvements achieved since Serber. Particle transport codes can then rely on such spallation models to 
treat the reactions between a light particle and an atomic nucleus with energies spanning from few tens 
of MeV up to some GeV. An overview of the spallation reactions modeling is presented in order to point 
out the incomparable contribution of models based on basic physics to numerous applications where such 
reactions occur. Validations or benchmarks, which are necessary steps in the improvement process, are 
also addressed, as well as the potential future domains of development. Spallation reactions modeling is 
a representative case of continuous studies aiming at understanding a reaction mechanism and which end 
up in a powerful tool. 
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1 Introduction 

The spallation reactions are one branch of the nuclear 
reactions. Usually they are defined as the interaction of 
a particle, most of the time a nucleon, with a nucleus. 
The projectile energy lies more or less between 100 MeV 
and 3000 MeV. During this interaction many nucleons are 
emitted and the residual nucleus can be very different from 
the target nucleus. 

Early studies of spallation reactions took place in 1947 
in the Berkeley cyclotron at the University of California 
[T]. Several experiments (items A8 and A9 of [2]) with 
beams of deuterons and alpha of energy 200 MeV and 400 
MeV, respectively, were able to demonstrate the multiple 
emission of neutrons and charged particles resulting in the 
production of a large spectrum of residual nuclei. Measure¬ 
ments of neutron spectra were also taken and compared, 
quite successfully, to the theory of R. Serber who was the 
first to describe, also in 1947, the process of the spallation 
reaction in two stages [3] . 


From a historical and semantics point of view, the spal¬ 
lation term seems to be first used by G. T. Seaborg in his 
thesis in 1937. The first publications mentioning the spal¬ 
lation term date back to 1948: O’Connor and Seaborg [4] 
(in the title) and Seaborg and Perlman ([5], page 587). 
For the record and as anecdote Serber didn’t mention the 
spallation term in his aforementioned article, but he titled 
it Nuclear Reactions at High Energies. The meaning of the 
word spallation is given below to avoid ambiguity later. 

It is usual to define the spallation reactions as nuclear 
reactions between a fast particle and an atomic nucleus. 
In recent decades physicists prefer to speak of protons of a 
few hundred MeV up to a few GeV on a nucleus, perhaps 
because during this time the spallation has been closely 
associated with the ADS (Accelerator-Driven Systems), 
for which the beam is made of protons. Actually the par¬ 
ticle type and the energy limits can be slightly different 
according to the researcher and/or to the period. More¬ 
over regarding the modeling we will see that the spallation 
reactions are a set of successive or competing nuclear pro- 
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cesses : binary collisions, pre-equilibrium, evaporation, fis¬ 
sion, binary decay, multifragmentation. In addition in all 
directions some extensions can be done, towards lower and 
higher projectile energies and using more complex projec¬ 
tiles. Thus since the goal is to describe nuclear reactions, 
like the spallation reactions, within the largest ranges of 
energy and of projectile type, the needed assumptions and 
approximations in the modeling must not only be kept in 
mind in order to know the validity domain, but they define 
this domain. In other words assumptions and approxima¬ 
tions define the reactions we are studying. Moreover, for 
sake of clarity, the definitions of the generic terms particle 
and energetic particle adopted in this paper, when they are 
used alone, are the following. Particle : lepton or hadron 
with a mass less than or equal to the a particle. There¬ 
fore all nuclei are heavier than the a particle. Energetic 
particle : particle having an energy such that it sees the 
atomic nucleus like a set of free nucleons. The lowest en¬ 
ergy of the incident particle is thus limited to around 150 
MeV . The upper limit is more difficult to determine, but 
it makes sense not to go far beyond the heaviest baryonic 
resonances implemented. 

The reasons and motivations for studying these reac¬ 
tions are manifold. Spallation reactions exist in several 
places and are linked to several domains. In space. Galac¬ 
tic and Solar Cosmic Rays (GCR and SCR) induce such 
reactions on spacecrafts and interstellar bodies like the 
meteorites. In the first case the knowledge of the reac¬ 
tions are mandatory for safety reasons. Radiation damage 
to the electronics and radioprotection of the astronauts 
must be under control in order to safely achieve the mis¬ 
sion [6]. In the second case the spallation reactions occur- 
ing in the meteorites can allow us to obtain information 
on their history via the residual isotopes produced inside 
provided we have sufficient knowledge of their production 
mechanism ([8] and references therein). In addition, GCR 
and SCR are partly responsible, via spallation reactions, 
of the abundance of Li, Be and B [7]. In accelerators, the 
spallation reactions are used most of the time either as 
neutron sources or as radioactive isotope beams. The term 
” spallation neutron source” is usually used to describe fa¬ 
cilities that provide neutron beams (which are eventually 
cooled to get the right wavelength) to probe matter. Such 
facilities have existed for a long time, for example SINQ 
(in Switzerland) [9] and ISIS (in UK) [10], and new ones 
have been recently built in the USA and in Japan, re¬ 
spectively SNS (operational in 2006) [TT] and JSNS (in 
2008) [12]. The next planned high-power neutron facility 
in Europe, the European Spallation Source (ESS) [13] will 
be built in Lund (Sweden) and should deliver first neu¬ 
trons around 2020, and more or less at the same time the 
China Spallation Neutron Source (CSNS) [M] should also 
be operational. However the neutron production via spal¬ 
lation reactions can also be used to control hybrid reac¬ 
tors. Such facilities, named ADS (Accelerator-Driven Sys¬ 
tem), consist of a nuclear reactor, with a spallation target 
surrounded by the fuel, and of an accelerator. The added 
spallation neutrons enable control of the reactor even with 
fuels containing actinides delivering only a low fraction of 


delayed neutrons. Since the reactor is subcritical by the 
prompt neutrons and the accelerator can be shut down 
quickly, such a facilty could help to transmutate nuclear 
wastes safely. The use of spallation neutrons in order to 
fission a nucleus can also be a way to deliver radioactive 
ion beams. A target, made generally of actinides and sur¬ 
rounding a spallation target, undergoes fission reactions 
via the spallation neutrons and the fission products, once 
extracted and selected, can be used as a secondary beam 
of exotic nuclei. A complementary way is to directly use 
the radioactive nuclei produced by spallation. According 
to the expected beam, one or the other of the methods has 
to be chosen as well as the material of the spallation target, 
to get the highest flux. This type of radioactive ion beam 
facility is called ISOL, for Isotope Separation On-Line, and 
we can mention three of them that are operating today: 
CERN-ISOLDE and GANIL-SPIRAL in Europe and TRI- 
UME in North America. In the more remote future a high 
intensity ISOL facility is planned in Europe: EURISOL. 
In the meantime, new facilities will be built and operated 
in the near future, for example HIE-ISOLDE (CERN), 
SPES (Legnaro) and SPIRAL2 (GANIL). Einally we can 
also cite another application, proton (or hadron) therapy, 
where the reaction destroys the cancerous cells, and the 
proton, unlike the usual radiotherapy, delivers only low- 
dose side-effects to surrounding tissue. 

In addition to those motivations driven by application 
of the reaction, the spallation reaction study is obviously 
interesting in itself. Understanding a reaction mechanism 
involving a particle, with an energy between 100 MeV and 
2-3 GeV, and any atomic nucleus (light and heavy), i.e. 
being capable of reproducing emitted particle spectra in 
energy and in angle, as well as the mass, charge and iso¬ 
topic distributions of the residual nuclei, is a real chal¬ 
lenge. This challenge in modeling is nevertheless neces¬ 
sary, because, due to the broad projectile energy range, 
to the large target nucleus spectrum and to the numer¬ 
ous output channels, each experimental data set can only 
be used to validate a part of the reaction mechanism. It 
is impossible to launch a program aiming at covering all 
types of spallation reactions. 

This paper aims at presenting the spallation reactions 
from modeling to application in parallel with the improve¬ 
ments over the last few years. It is divided into four parts: 
modeling, validation, application and future development. 
Most of the time publications dealing with spallation re¬ 
actions address only one of the first three parts, and when 
two of them are addressed, validation is included, but with 
a minor and ancillary role. The first section, modeling, is 
focussed on the outlines. A review paper could be written 
only on spallation reaction modeling. We then aim at cov¬ 
ering only the main aspect of the modeling from the basic 
ideas up to the computational codes. Details are some¬ 
times given to show either the various domains involved 
in the modeling or the difficulties met and the assump¬ 
tions needed. Eor example two models will be explained 
more in depth, INCL4 for the first stage and Abla for 
the second stage of the reaction, because, first, these two 
codes have been extensively used in the part dedicated to 
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applications and, second, they are good examples of mod¬ 
eling evolution. Nevertheless, as far as possible, numerous 
other spallation codes are discussed in the second section 
devoted to validation. Whether combined to modeling or 
to application, validations are usually addressed in a spe¬ 
cific purpose. Here validation is discussed in a more gen¬ 
eral way with a brief history of various benchmarks carried 
out in the two last decades, with an overview of the bench¬ 
marks done by the transport code developers and with 
the results obtained within the Benchmark of Spallation 
Models, an IAEA project where the best known spallation 
models were involved, comparing their calculation results 
to an experimental data set chosen in order to cover all 
aspects of the spallation reaction, in input channels and 
output channels. The third section gives four examples 
of the use of spallation reaction modeling: the EURISOL 
project, the MEGAPIE target, the ESS facility and the 
cosmogenic nuclide production in meteorites. The goal is 
to discuss respectively designing and optimisation, pre¬ 
dictive power and reliability, feasibility and uncertainty 
estimate, and the use of modeling to mitigate lack of ex¬ 
perimental data. Einally the fourth section is devoted to 
the topics less studied up to now, like e.g. the pion and 
strangeness production, especially towards the high ener¬ 
gies. 

This paper does not aim at being complete, giving all 
details in all topics. Such an objective is out the scope of 
a review article, and that has been successfully addressed 
in the book entitled Handbook of Spallation Research [15], 
written by Eilges and Goldenbaum. The goal of this pa¬ 
per is rather to point out the reasons for the successes 
of the simulations without hiding the remaining deficien¬ 
cies and to give references to help the reader go to into 
further detail. This is also the opportunity to show the re¬ 
sults achieved by physicists involved in spallation reaction 
study over several decades. 

2 Modeling and associated codes 

The basics of modeling spallation reactions have been stated 
by R. Serber in 1947 [3]. Unlike the low-energy reactions 
where there is formation of a compound nucleus, which 
then de-excites, the high-energy reactions, say above 100 
A.MeV, cannot form directly in general a compound nu¬ 
cleus due to the large mean free path of the incident nu¬ 
cleon. This large mean free path, thus associated with a 
small de Broglie wavelength, allows the incident nucleon to 
see the nucleons in the nucleus as if they were free and to 
start a cascade of binary collisions. Different scenarios for 
this cascade are possible depending on the energy and the 
impact parameter of the incident particle going from ejec¬ 
tion of a single nucleon, taking with it all of the incident 
energy, to the capture of the projectile leaving the nucleus 
in a state of strong excitation. Once excited, the nucleus 
enters a second and slower phase, the de-excitation. Here 
again different scenarios compete according to the mass, 
excitation energy and angular momentum of the remnant 
nucleus. The first and rapid phase is of about 10“^^ s 
and the second in the order of 10“^^ s. In addition to 


these two phases sometimes included is a third one named 
pre-equilibrium. This step is actually an intermediate step 
since it deals with the transition between cascade and de¬ 
excitation and more precisely how the cascade is stopped. 
The need of this additional phase is then strongly con¬ 
nected with the cascade modeling. 

The codes associated with the modeling are of course 
the codes handling each phase of the spallation reactions, 
but there exists also others associated codes: the particle 
transport codes. These codes aim at simulating interac¬ 
tions of particles with a massive target, i.e. which has a 
spatial extension and therefore is the seat of successive re¬ 
actions due to the secondary particles. They usually use li¬ 
braries to handle the transport of particles, but, according 
to the type of the particle and its energy, reaction models 
are needed when the experimental or evaluated data are 
unavailable. This is typically the case of spallation reac¬ 
tions. Although some libraries have emerged [16][T7|, they 
have been built thanks to models. 

This section is divided into three parts. The first two 
ones deal with the two above mentioned phases of the 
spallation reaction (the pre-equilibrium being addressed 
with the cascade). We will give the main aspects, the basic 
ideas, and details only when necessary to understand the 
topics discussed in the following sections of this paper. 
The last subsection is devoted to transport codes related 
to the spallation regime, because these codes are most of 
the time the link between modeling and applications. The 
transport codes need models, but in return the models had 
to be upgraded to be used optimally in transport codes. 


2.1 Intranuclear cascade (INC) - First phase 

Modeling of the first phase consists of simulating the time 
evolution of a system of nucleons in a nucleus possibly 
in interaction with each other. Several articles have been 
published on the basic foundations of spallation reactions 
and include among others Aichelin [18] , Bunakov [19] and 
Cugnon [20]. Aichelin and Bunakov demonstrate in two 
different ways how the problem moves from quantum- 
mechanical equations to transport equations. Bunakov in 
addition explains that the intranuclear cascade (ING) mod¬ 
eling is nothing else than a Monte-Carlo method to solve 
the transport equations. This point is taken up in m 
by Cugnon to explain the good results obtained by some 
codes using this modeling. In the following we summarize 
the very main steps to get the master equations. In this 
way we will show what are the difficulties to be solved 
and the approximations to be made to be able to simulate 
spallation reactions. 

The time evolution of a quantum system is given by 
the Schrodinger equation, using the wave function, or by 
the Liouville-Von Neumann equation (eq. [^, using the 
density matrix p, 

tn^ = [H,p], ( 1 ) 

where H is the Hamiltonian. 
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The goal being to solve this equation, a first step is to 
use reduced density matrices (eq. , which are the density 
matrices of a subsystem of k nucleons, 


Pi,...,k = 


A\ 

{A-k)\ 




( 2 ) 


such that = P, with A the total number of nu¬ 

cleons. Doing so the Liouville-Von Neumann equation is 
replaced by a system of A coupled equations (eq. [^, so- 
called BBGK\0 hierarchy. 


..■*> 

^ k k 

N [^i, p{t)\ + 2 ^ [rij, /?(i)] + 

I *=1 i,j 

\ 

k 

i=l 

( 3 ) 



Such a system does not simplify the solution of the A- 
nucleons system, however it allows to reduce the chain 
through physical approximations and so make the system 
solvable. In our case collisions between two particles have 
to be taken into account and they occur in a nucleus that 
can be characterized by its potential. Thus only the first 
two equations (eqs. and have to be considered, with 
the two interacting particles (1 and 2) and a third partner 
(3), which represents the nuclear potential, 

= [Ti, Pi{t)] + ^T^2) [Vi2, Pi,2{t)], ( 4 ) 

= [Ti + r2i Pi, 2 ( 0 ] + 

II (5> 

2 [^12,Pl,2(^)] + [^13 + k23,Pl,2,3(^)] • 

Ti and Vjk are respectively the kinetic energy and inter¬ 
acting potential operators of the particles i, j and k. 

Through the use of the Wigner transformation, /, of 
the density matrices, 

f{r,p) = J(fr' p{r +'^,r - 

and the use of three main assumptions, which are (as men¬ 
tioned in [2Q]h 

• the two-body correlations are unimportant in the col¬ 
lision term, 

• the mean potential is a smooth function on the spatial 
extension of the particle wave functions, and 

• the collisions are independent, that is the scattering 
wave function becomes asymptotic before the next col¬ 
lision takes place, 

^ BBGKY: Bogolioubov-Born-Green-Kirkwood-Yvon 


the transport equation, also called BUU or VUlj[^ equation 
[ 21 ], can be written as 


^ + V.Vrfl - VrU.Vph 


-J 


(fp2d^Pzd^PA 


s^ipl +P 2 -P 3 


(27r)6 
PA)S{ei + €2 - 


avi2 


es 


£ 4 ) 


[/ 1/2 (1 - /a) (1 - Ja) - fsfA (1 - /i) (1 - / 2 )] • 


fi are the distribution functions (1 and 2 for the two collid¬ 
ing particles, 3 and 4 for the two outgoing particles), avi 2 
the collision term (most of the time based on experimen¬ 
tal data), the fifj and {1 — fi) terms are respectively the 
overlapping of input/output distributions and the Pauli 
principle, and finally the delta functions are the momen¬ 
tum and energy conservations. 

Doing so Bunakov m not only explained how it was 
possible to move the problem from quantum-mechanical 
equations to transport equations, through several assump¬ 
tions matching spallation reactions, but also demonstrated 
the basic foundations of the intranuclear cascade (INC) 
modeling. 

From the practical point of view the different steps 
encountered during the intranuclear cascade in all mod¬ 
eling codes are: initialization, transport, event treatment 
and stopping. Since all these details are explained for each 
available model in articles and/or user’s manuals (a list of 
references is given in section]^, only the main features are 
described below. Nevertheless three papers must be cited 
here. The two articles of Metropolis et al [22l |23] , where 
the first calculation results about spallation reactions us¬ 
ing a compute]0 were published, and the one of Bertini 
[25] . where an improved version of the previous model is 
used. Those works are the basis of modern codes. 

Both projectile and nucleus have to be initialized. Re¬ 
garding the projectile, type and energy are known, but its 
impact parameter is taken randomly and coulomb devia¬ 
tion considered. If the projectile is a composite particle, its 
structure must be given in the same way as for the target. 
The target nucleus is defined by its mass, its charge, the 
potentials felt by the particles, the momentum of each nu¬ 
cleon (most of the time a Fermi gas distribution is used), 
and the spatial distribution of the nucleons. Two ways ex¬ 
ist to define this distribution. Either the distribution is 
continuous, often several concentrical density regions, or 
discrete, i.e. positions are sampled in a Wood-Saxon dis¬ 
tribution, for example. 

Two different methods are applied to move and follow 
the particles participating to the cascade. With the time- 
like transport all the particles are followed at the same 
time, while with the space-like transport the patrticles 
are followed ones after the others. 

Three events exist during the intranuclear cascade: 
collision, resonance decay and refiection/transmission at 
the nucleus surface. Collisions can be elastic or inelas¬ 
tic. Most of the time experimental data (cross-sections) 

^ BUU: Boltzmann-tihling-Uhlenbeck; VUU: Vlasov-Uhling- 
Uhlenbeck 

^ The first calculations by hand were done by M. L. Gol- 
berger in m- 
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are used to define interaction probabilities first and sec¬ 
ond what are the output products and their characteris¬ 
tics (types, energies, momenta), each selection being done 
randomly. When necessary Pauli blocking is taken into ac¬ 
count. Cross-section parameterizations, number and type 
of collisions and of resonances taken into account, and the 
way to apply Pauli blocking are different from one model 
to another; this can already explained the broad spec¬ 
trum of results that will be shown in section Decay of 
the resonances is not as simple as it seems. It is not only 
the production of two particles (Pauli blocking included), 
but also the right behavior of the resonance. This parti¬ 
cle has a lifetime that can be more complicated than the 
inverse of the decay width. For example, a phase-space 
factor has been added to the Z\(1232) lifetime, because 
this resonance is a irN system at low energy. More details 
are given in [26], where the Z\(1232) has been studied. 
Moreover, experimental data for the other resonances are 
sometimes scarce, and, due to their widths and proximity 
(in energy), some can overlap and make difficult a clear 
treatment. This explains, in addition to the energy do¬ 
main, why sometimes only the Z\(1232) is implemented. 
The third event happens when the particle reaches the 
surface; it can be emitted or reflected. To be emitted the 
particle must be energetic enough, i.e. be able to over¬ 
come nuclear and coulomb potential. Experimental results 
on composite particle spectra display high energy tails, 
which indicate emission of composite particles during the 
intranuclear cascade phase. However the basis of INC are 
clearly the treatment of the transport of nucleons with 
their two-body interactions, i.e. without clusterization. So, 
up to now, the only way to produce composite particles 
during the cascade is to add a coalescence model. Before 
leaving the nucleus, a nucleon can drag and aggregate one 
or more nucleons, close enough to it in space and momen¬ 
tum. This procedure extends the INC applicability to a 
satisfactory level. 

Finally different criteria are used to stop the cascade 
and to start the de-excitation phase of the remnant nu¬ 
cleus. We can mention three of them: cutoff energy, stop¬ 
ping time and deviation from an optical absorptive poten¬ 
tial. The first one is the simplest. When its kinetic energy 
falls below a given cutoff energy (depending on the type 
of the particle or not) the particle no longer contributes 
to the cascade and, once there are no more participants, 
the cascade is stopped. The stopping time criterion, de¬ 
fined and used in INCL4 m, is based on physical fea¬ 
tures. The two stages, cascade and de-excitation, are two 
different mechanisms (duration, characteristics of the ini¬ 
tial and final products) and the hot remnant nucleus, de¬ 
exciting during the second phase, can be seen as a com¬ 
pound nucleus, i.e. with no memory. In [27] the authors 
of INCL4 describe the result of their study that ends up 
in an equation giving this stopping time versus the tar¬ 
get mass. Figj^ shows the time evolutions of four observ¬ 
ables used to determine this time. The last and third cri¬ 
terion, used in CEM03 [29], is based on optical potentials. 
The cascade is stopped when the deviation between the 
imaginary parts of the calculated optical potential (from 


interaction cross-sections and Pauli blocking effects) and 
of a phenomenological optical potential exceeds a chosen 
value. The cascade formalism is then considered inappro¬ 
priate. This criterion is described in [29] . In this paper the 
authors explain the difficulties to find the right value to 
assess this deviation, which governs the switch to a pre¬ 
equilibrium model. 

Such a pre-equilibrium model is sometimes used be¬ 
tween the cascade and the de-excitation phases. Several 
versions exist, but almost all are based on the exciton 
model developed by Griffin [30]. According to the use or 
not of this intermediate phase, the duration of the cas¬ 
cade is obviously different or, maybe more correctly, mass, 
charge and excitation of the remnant nucleus are larger, 
if this phase is called. While some intranuclear cascade 
models need such pre-equilibrium models to improve their 
capability, this is not the case of some others. This will be 
shown in section and this topic is related to the low 
energy interaction during the cascade. 

After mentioning in [20] the three main assumptions 
needed to draw the master transport equation, J. Cugnon 
points out also that the INC models do more than solving 
this equation. In particular the first assumption, the two- 
body correlations are unimportant in the collision term, 
is not made, since the previous collisions influence the 
next ones. Thus collision after collision n-body correla¬ 
tions are propagated leading to reliable prediction of the 
spectral characteristics of secondary particles. Bunakov al¬ 
ready discussed the domain of validity of INC model from 
its basic assumptions. He considered that INC method 
can be applied with a target mass above A^IO and with 
an incident projectile energy above about ten MeV, and 
added that from the point of view of its basic principles 
the INC model is not inferior to the other reaction models, 
say Glauber^s model, the optical model or DWBA. Never¬ 
theless the INC method has obviously its own limits, such 
as the need of a smooth function for the potential that 
makes the use of INC difficult to treat Nucleus-Nucleus 
reactions and some justifications are at least necessary. 
Still from the basic assumptions, Bunakov suggests some 
modification to extend or improve the INC capabilities. 
We can cite the treatment of the nucleon mean free path, 
which should be calculated carefully for the low energy 
particles (< 100 MeV), and the use of DWBA calculations 
for the first collisions to improve the threshold behavior. 
Since then INC models have been refined and their reli¬ 
ability significantly improved, thus we mention hereafter 
two studies related to Bunakov suggestions. 

In their paper [31], J. Cugnon and P. Henrotte in¬ 
vestigate the low-energy limit validity of the intranuclear 
cascade on a purely empirical basis. They compare their 
INCL4.3 model to experimental data with projectile en¬ 
ergies from ~25 MeV to ^250 MeV and sometimes to 
other models dedicated to this energy domain. The first 
observation was the rather good results well below the of¬ 
ten quoted lowest limit and the second was that INCL4.3 
could compete with more dedicated models, even if there 
was room for improvement. They then study two main 
features: the Pauli blocking implementation and the num- 
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Fig. 1: Time evolution of four observables during the in¬ 
tranuclear cascade INCL4. Nucleus excitation energy is 
plotted at the top-left, mean kinetic energy of the emitted 
particles at the top-right, time derivative of the excitation 
energy at the bottom-left and momentum asymmetry of 
the participants at the bottom-right. The arrows point out 
the time that can be considered as the end of the cascade. 
Reprinted figure with permission from A. Boudard [27]. 
Copyright (2014) by the American Physical Society. 


ber of collisions. Pauli blocking can be accounted for via 
different ways. In this study they used a ’’strict” and a 
” statistical” version. ” Strict” means that collisions are al¬ 
lowed only if the momenta of the outgoing particles are 
larger than the Fermi momentum. ” Statistical” version is 
based on the phase-space occupation probabilities around 
the outgoing particles and a comparison with a random 
number. The conclusion was the best results were obtained 
using both of them, the ’’strict” Pauli blocking for the first 
collision and the ’’statistical” one for the others. Accord¬ 
ing to the way positions and momenta are assigned to the 
nucleons to prepare the target nucleus, and taking into 
account the population evolution, the best implementa¬ 
tion was a priori not so easy to define. Nevertheless it is 
clear that the Pauli blocking plays a significant role in the 
mean free path. Regarding the number of collisions, it is 


defined as the level of the collision in the cascade, in other 
words as the incrementation of the number of collisions 
undergone previously. This study combined with a com¬ 
parison of other models, most of them using exciton states, 
shows first that INCL gave fairly good results compared 
to these more elaborate models in this low-energy domain, 
but also that the collision number was closely related to 
the exciton number. In addition to the comparisons with 
experimental data, this latter observation is a clue that a 
pre-equilibrium phase is not mandatory in the spallation 
modeling, provided the INC model treat carefully the low 
energy collision. 





Emitted neutron energy [MeV] 



Fig. 2: Double differential cross-sections and thick target 
yields (respectively top and bottom panel) at =3.5° in 
^Be(d, xn) for 18 MeV. Experimental data [33] are repre¬ 
sented by closed circles, and results calculated by Parti¬ 
cle and Heavy Ion Transport code System (PHITS) with 
JQMD, INCL, and INCL with DWBA are denoted by dot¬ 
ted, dashed, and solid lines, respectively. Figures drawn 
from [32]. Copyright @ 2014 Elsevier B.V. 


The use of a DWBA model to improve behavior at 
the low-energy limits has been investigated recently. This 
study combines a version of INCL, INCL4.5, to a DWBA 
method [32], because INCL treats well the incoherent pro¬ 
cesses and DWBA the coherent ones. Details on the two 
models can be found in the references, but what has to 
be mentioned here is the need for fitted parameters in 
DWBA, especially for the optical potential used. Thus 
this study focused on (d^xn) spectra (and also {d^xp)) 
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on ^Be, and targets at incident energies rang¬ 

ing from 10 to 40 MeV. From the technical point of view 
the use of INCL or DWBA is chosen randomly according 
to the ratio of the reaction cross-sections. We plot 

results in Fig.|^ one is the elementary production of neu¬ 
trons at a forward angle (upper part) and the second one 
a neutron thick target (2.5 mm) yield (lower part), INCL 
and DWBA being implemented in the PHITS transport 
code. We clearly see the two domains: INCL4.5 reproduces 
well the neutrons at half the incident energy of a deuteron 
and DWBA the high energy part (the first colision having 
then a heavy weight). 

2.2 De-excitation - Second phase 

As previously stated, spallation reactions are divided in 
two mechanisms: the intranuclear cascade and the de¬ 
excitation. This can be illustrated with Fig. where neu- 



Fig. 3: Neutron double differential cross-sections for the 
reaction p(1200 MeV) + W. Reprinted figure with per¬ 
mission from S. Leray [ml where details can be read. 
Copyright (2014) by the American Physical Society. 

tron spectra are plotted for several emission angles. The 
high-energy part (cascade) is strongly focused in the for¬ 
ward direction and characterized by two peaks, the elastic 
one (sharp) and the quasi-elastic one (broader) due to 
the Z\(1232) resonance, whereas the low energy part (de¬ 
excitation), much more important, is isotropic. The large 
amount of emitted nucleons during the de-excitation stage 
can also be seen in mass and charge distributions of the 
residues, as shown in Fig.|^for the reaction ^^Fe+p, where 
residue mass distributions measured at five projectile en¬ 
ergies, from 300 A.MeV up to 1500 A.MeV (here inverse 


kinematics has been used), are plotted. The upper part 
of Fig. shows also that the production of the light nu¬ 
clei increases with the projectile energy, in such a way 
that, as seen on the lower part, the yield of some very 
light residual nuclei becomes significant. This latter ob¬ 
servation already indicates that the emission of nucleons 
is not the only way to de-excite. Another complexity of 



A 


I ' 1 ’ I ^ 

100 - 



Fig. 4: Mass distribution of the residual nuclei in the 
spallation reaction ^^Fe+p. Above at five different beam 
energies (Reprinted figure with permission from C. Vil- 
lagrasa nHi ) and below at 1000 A.MeV (Reprinted fig¬ 
ure with permission from P. Napolitani [182]). Copyright 
(2014) by the American Physical Society. 


the de-excitation phase arises when the target mass in¬ 
creases as shown with Fig.[^ On this charge distribution, 
in addition to the decreasing cross-section with decreas¬ 
ing residual charge, on the high Z part, a bump centered 
around the half of the target charge appears. This bump 
indicates a high-energy fission process. 

Once the nucleus target reaches an equilibrium state, 
the intranuclear cascade vanishes and is followed by the 
de-excitation of this compound nucleus (remnant) char¬ 
acterized by its mass, charge, excitation energy and an¬ 
gular momentum. The most common mechanisms used to 
describe this phase are the particle evaporation and the 
fission process, but others have been proposed either as al¬ 
ternative ways or as additional channels to cover the whole 
mass (and charge) distributions. Thus we report hereafter 
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Usually the state density is defined with the Fermi gas 
model using a canonical ensemble rather than a micro- 
canonical one to simply the formalism. Then a new param¬ 
eter arises, the level density parameter a, which relates the 
excitation energy and the temperature: E"* = aT^. Thus, 
considering the nucleons distributed on equidistant levels 
at the Fermi energy (place of the excitation energy), the 
state density, with no angular momentum dependence, can 
be written 


p{E*) 


^2VaE* 


12 


(7) 


and with angular momentum dependence J 


Fig. 5: Charge distribution of the residual nuclei in the 
spallation reaction ^^^Pb+p at 1000 A.MeV (Reprinted 
figure with permission from A. Boudard EH). Copyright 
(2014) by the American Physical Society, Data |37] are 
the black points. 


the five mechanisms possibly encountered in the spallation 
modeling for the de-excitation step: evaporation, fission, 
binary breakup, multifragmentation and Fermi breakup. 
The same approach is used in all these processes, that is 
the statistical approach. The number of microstates of a 
given macrostate (here the remnant nucleus) determines 
the probability of a channel and the final channel is chosen 
randomly. 


2.2.1 Evaporation 

Understanding the evaporation of particles means to be 
able to determine: i) the type of the emitted particle or 
more precisely the ratios between all types , ii) how many 
of these particles are emitted for each type, and iii) the 
energy spectra of the emitted particles. To answer these 
questions the emission probability of each type of particles 
as a function of the kinetic energy must be calculated. This 
is done starting with the following assumptions: 

• the nucleus is sufficiently excited to apply a statistical 
description, and 

• the emission rate is equal to the rate of absorption 
(detailed balance principle). 

Doing so one obtains the emission probability per units of 
time and energy of a particle {part) with a kinetic energy 


part 


(^) = 


p^{E:) 


{2s + l){47rp^ dp/hE a c{s), (6) 


p(E*,J) = (2J + l)(^ 


3/2 


^/d e 


2vEU 


12 m 


( 8 ) 


where I is the nucleus inertia moment and U = —Erot^ 
with Erot Ihe rotational energy. The point is then the ex¬ 
pression of the level density parameter a. Numerous works 
deal with a and it is impossible to draw here a list and 
even a summary of the different values and/or formulae 
used. Nevertheless the knowledge of this parameter has 
been improved since the first estimates. In [38] WeiBkopf 
gives values for a, in [39] Gilbert and Cameron take into 
account the pairing and the shell effects, and later these 
effects have been refined in [40] and [41]. The deforma¬ 
tion of the nucleus and the diffusivity of its surface have 
been investigated in [424144] . even if the temperature de¬ 
pendence can be neglected in the spallation domain, this 
has been studied in m and the collective effects are in¬ 
cluded in the level density parameter used in the Abla 
code [46] . 

The second basic ingredient for determining the width 
or probability of emission of a particle is the capture cross- 
section as a function of its kinetic energy 5. Various ex¬ 
pressions can be found in literature, from the simplest one 
to the most sophisticated. The simplest formula is the one 
used in 1937 by WeiBkopf [38], but taking into account the 
Coulomb barrier Bcoui of charged particles 


(Jc{^) = 'kR^{1 — 


B. 


coul 


(9) 


with R, radius of the nucleus. In 1940 WeiBkopf and Ew¬ 
ing improve the formula in m- It takes into account the 
angular momentum by a decomposition into partial waves 
and the cross-section becomes 






( 10 ) 


where k is the wave number and I the angular momentum, 
where/9^(E^ )(/9 j(Ej)) is the state density of the initial (final) WeiBkopf-Ewing split their Qi{£) factor into two parts, 
nucleus with E* {E}) the excitation energy of the initial(final)The first one is related to the incoming wave, P; (e), which 


nucleus, (Anp^dp/h^) the phase-space of the emitted par- 
tide with a momentum p, s its spin, and (Jc{s) the capture 
cross-section of the particle by the final nucleus. The two 
main ingredients are the state density and the capture 
cross-section. Their choices largely explain the differences 
between the codes. 


is the penetration probability into the surface of the nu¬ 
cleus, and the second one related to the energy exchange 
probability between the particle and the nucleus, ^(e), 
called sticking probability. Pi{s) contains the barrier for 
charged particles (actually the Coulomb barrier and the 
centrifugal barrier) and, regarding ^i{e), only the average 
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value on I was available. A very similar formula is given 


by Hauser-Eeshbach [48]. 



(11) 

The Qi (e) factor is replaced by the penetrability or trans¬ 
mission factor Ti{e), which has become more sophisti¬ 
cated. Eor example it can be written as 


(12) 

with 

(13) 


where /i is the reduced mass of the system and V{I, r) the 
potential taking into account the Coulomb potential, the 
centrifugal potential and the nuclear potential. While the 
centrifugal potential is clearly known and Coulomb poten¬ 
tial also, although refinements exist, the nuclear potential 
is more complicated and different expressions exist such as 
those of Huizenga [49] or Bass m- Whatever the formula 
used for the capture cross-section, it amounts mainly to 
calculating the V{l^r) potential. 

When using experimental data to fit the formula of 
the inverse or capture cross-section, one has to pay atten¬ 
tion to the fact that the compound nucleus is not in its 
ground state, but heated. This is generally not a problem, 
except at low energy. Some other aspects must be taken 
into account, like the tunneling effect. This point is obvi¬ 
ously important for the low-energy light charged particles, 
as explained in the section describing AblaOT in m- 

The formalism of evaporation is well posed, especially 
for the neutron, but the expressions of its ingredients may 
vary and explain the differences in results according to the 
models, mostly for charged particles. For these latter par¬ 
ticles it is worth mentioning that evaporation models often 
also emit light nuclei. While typically the same formalism 
is used for nucleons up to alpha, the strategy can change 
when the number of channels increases (emission of any 
type of nuclei), because the number of emission widths 
or probabilities increases significantly and the computa¬ 
tion time also. The need to evaporate nuclei with a charge 
higher or equal than 3 has been seen in Fig. [^ This is 
however not the only way to generate such nuclei. Other 
mechanisms, such as binary decay or multifragmentation, 
are also used to do the job. 

Finally the gamma evaporation is a special case. Gam¬ 
mas do not compete with the other channels, except at the 
end of the de-excitation stage, when the nucleus is much 
colder. Then neutrons and gammas compete. Most of the 
codes include the statistical emission of gammas through 
the giant dipole resonance. However this treatment ends 
up not always with the splitting into the ground state 
and the metastable states. Nevertheless we can mention 
a recent model, using in particular the Evaluated Nuclear 
Structure Data File [52|, that has been developped within 
the code PHITS and gives interesting results on isomer 
production [53]. A comparison between the codes provid¬ 
ing the isomer production should be interesting, because 
this domain has not been studied very much up to now. 


2.2.2 Fission 

When excited a heavy nucleus can also undergo fission. 
This process competes with all evaporation channels and 
even with all de-excitation channels. The treatment of the 
fission process is indeed divided in two separate parts: cal¬ 
culation of the fission width (or probability), which gov¬ 
erns the competition, and, if fission occurs, determination 
of the characteristics of the two fission products (mass, 
charge and excitation energy). 



Fig. 6: Dependence of the barrier heights on /A. The 
points are extracted from the analysis of nuclear fissility 
done in ESI- The curve is liquid drop model calculation re¬ 
sults. Figure from [55] (Courtesy of V. Muccifora). Copy¬ 
right © 2014 Elsevier B.V. 

Two parameters explain the instability of excited nu¬ 
clei against fission: the fissility parameter and the 

angular momentum J. The fissility parameter describes 
the competition between the repulsive Coulomb energy 
(oc Z ‘^!and the attractive surface energy (oc ^2/3) 
(see EigJ^. The higher the fissility, the lower the fission 
barrier. The angular momentum increases the deforma¬ 
tion of the nucleus and then decreases the fission barrier. 
The results obtained in the work of A. Sierk on a macro¬ 
scopic model for rotating nuclei [54] are used by most of 
the fission models in the spallation domain. 

Eission, unlike evaporation, is a collective phenomenon. 
So passing from a state of individual degrees of freedom 
to a state of collective degrees of freedom takes time, the 
transition time to reach the highest fission probability. The 
random movement of particles in a dissipative medium is 
managed by the Langevin or Eokker- Planck equations. If 
the stationary equations are solvable analytically, this is 
not the case for the dynamical ones and long computation 
times are needed. As is often done to save time, approxi¬ 
mations are required, but the results are convincing [56] , 
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Two consequences of what has been stated previously: tak¬ 
ing into account the viscosity via the stationary equations, 
the fission width is lower than expected, and second the 
maximum fission width value needs time to be reached 
(dynamical effects). The width can then be expressed as: 




where 

• ^Bw is the Bohr-Wheeler fission width, 

^ 271 pi{E*,J) J/^-ddieiE* - Bf{J) - e)de, (14) 

with Psaddie{E* — Bf{J) — e) the level density at the 
saddle-point and Bf{J) the fission barrier, 

• K the Kramers factor (= 1 + “ (^) ’ 

(3 is the reduced dissipation coefficient and cjo the fre¬ 
quency of the harmonic oscillator describing the po¬ 
tential at the saddle-point deformation) , and 

• fp(t) a function simulating the dynamical dissipation 
effects. 

Once the fission channel is chosen, characteristics of 
the fission fragments must be determined. In a few words 
we can say that experimental data are usually used for 
the determination of masses and charges, often assuming 
Gaussian distributions, and taking care to preserve the 
initial mass and charge. The excitation energy is shared 
in respect to the mass ratios, and, if the kinetic energy of 
each fragment is based on their Coulomb repulsion, ad¬ 
justments on experimental data is often needed also. One 
of the most reliable models developed on this topic is the 
one of J. Benlliure m- Once again, all models used to 
model a part of the spallation process must be fast enough 
to give results in a reasonable time, thus, combining this 
constraint with the goal of building a predictive model, 
J. Benlliure worked out a semi-empirical nuclear-fission 
model. 

The yield, Y{E^N)^ of a given fission fragment, here 
characterized by its number of neutrons, is the ratio of 
its number of transition states above the potential, at the 
fission barrier, to the sum of these numbers for all possible 
fragments: 


Y{E,N) 


p.(U)dU 


p is the density of state, E the excitation energy, V(N) the 
potential energy at the fission barrier and the number 
of neutrons in the fissioning nucleus. It turns out that this 
yield can also be written in an approximate way. 


Y{E,N) - 


with a a level density parameter simply taken equal to A/8 
in this case and £’* = E — V{N) the available energy above 
the potential. Therefore determining the distribution of 


fission fragments amounts to define the potential V{I^. 
The formula used is the following (colors related to Fig.I^: 

V{N) =Csym{N^d2 - Nf+ 

5Ui+ (15) 

Casym(i){Ni-Nf 

and an example is given for in Fig. 0to illustrate the 
different contributions. 



Fig. 7: Potential energy at the fission barrier for 
versus the neutron number of one of the fragments. Colors 
of arrows match the contributions in equation The 
original figure comes from [57] (Courtesy of J. Benlliure). 
Copyright © 2014 Elsevier B.V. 


The macroscopic part of the potential is represented by 
a symmetrical parabola with the Csym parameter and the 
microscopic part by the shell effects characterized by the 
two positions Ni {i = 1,2), the respective strengths SUi 
and the two curvatures Casym{i) • These parameters are ad¬ 
justed on experimental data. Finally a random sample is 
used to determine the N of the fragment. Other character¬ 
istics (charge, excitation energy) are obtained by including 
considerations of conservation, like the ratio of the num¬ 
bers of neutrons and protons in the fissioning nucleus, but 
may take into account some effects seen experimentally as 
the polarization charge in the case of the asymmetric fis¬ 
sion. To reproduce at best the experimental facts, features 
are defined at the most appropriate time. Thus the mass 
asymmetry is defined at the saddle point, low variation 
between the saddle point and the scission point, while the 
N/Z ratio is defined at the scission point, because it can 
vary even after the saddle point. 

Most of the other fission models used in spallation 
modeling are based on the Atchison’s fission model [58] 
with refinements and new parameterizations, as the model 
implemented in PHITS by Furihata [59] . 


2.2.3 Binary breakup 

Evaporation of particles and high energy fission can ex¬ 
plain almost all the charge and mass distributions, except 
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for intermediate mass nuclei, as previously mentioned with 
the bottom part of Figj^ Several means can be used to 
complete this description, like the extension of evaporation 
to heaviest nuclei than a, taking into account their excita¬ 
tion energy, or the multifragmentation process addressed 
in the next subsection. Another idea for the creation of 
these nuclei is the extension of the binary breakup. Evap¬ 
oration can be seen as a highly asymmetric binary decay 
and high energy fission as a symmetric binary decay. For 
Moretto [60], [...] there is no need to eonsider the two ex¬ 
tremes of this distribution as two independent proeesses. 
Rather, one would eonelude, fission and evaporation are 
the two, partieularly (hut aeeidentally) obvious extremes 
of a single statistieal deeay proeess, the eonneetion being 
provided in a very natural way by the mass asymmetry eo- 
ordinate. He then generalized the model of Bohr-Wheeler 
fission by adding an extra dimension: the mass asymme¬ 
try. The decay width of the extended binary breakup is 
then 

r,iE-.j)d,= _ 

where y is the mass asymmetry, py its conjugate momen¬ 
tum and niy the associated inertia. A similar, but simpler, 
formalism is used in the de-excitation code GEMINI, de¬ 
scribed in m- Charity applied a simplified formula for 
the width where the mass asymmetry is contained in the 
barrier 

rzAE\J)= 

2 ^- Psaddle{E* - BzAJ) " ^ ^ 

where Bz^a{J) = Ba{J) + -f AEcoui — — 6P, with 

AM and AEcoui the mass and Coulomb energy correc¬ 
tions due to the mass asymmetry, 6W and 6P respectively 
the shell and pairing effects. 

In practice, this binary decay process, useful in some 
codes to produce the intermediate mass fragments, are 
still competing with the evaporation and fission channels, 
whose the well established formalisms give good results. 


2.2.4 Multifragmentation 

We previously mentioned that intermediate mass isotopes 
appear when the excitation energy increases. Their pro¬ 
duction can be taken into account by an extensive evap¬ 
oration, i.e. beyond the a particle or by a generalized bi¬ 
nary breakup. Another physical process is suggested, the 
mult ifr agment at ion m, i e. the breakup of the target nu¬ 
cleus in several nuclei prior to its thermal equilibrium. 
The idea underlying this formalism comes from the time 
equilibration of the excited nucleus between two succes¬ 
sive emissions. This relaxation time becomes compara¬ 
ble to the time between two emissions of particle when 
excitation is greater than about 3 MeV, and thus the 
models based on a thermal equilibrium nucleus would no 
longer be justified. This relaxation time is equal to 2R(cs, 
where R is the radius of the nucleus and Cs the speed of 


sound in the nucleus. Estimates of the latter parameter 
can be found in reference [62]. To treat this multifrag¬ 
mentation different formalisms have been developed and 
we can cite, according to a ranking of J. P. Bondorf [63] . 
TDHF0 and qmeQ models for the microscopic dynamics 
and BUU/VUU for kinetic models. Another formalism is 
preferred in the case of spallation reactions with high ex¬ 
citation energies, i.e. with numerous degrees of freedom: 
the statistical approach. 

Briefly, the statistical model seeks to give a weight to 
each partition of the system and then one is randomly 
selected. Microcanonical, canonical and grand canonical 
ensembles can be used according to the knowledge of the 
system, to the possible approximations and to the neces¬ 
sary time to achieve reliable results. Thus, taking the case 
of the SMMj model, described in m, the number of par¬ 
titions quickly becomes huge (for 100 nucleons, ~ 2.10^, 
and for 200, ^ 4.10^^ jSl]), and simplifications based on 
physical findings are made. For rather low excitation ener¬ 
gies the microcanonical ensemble is used, but considering 
only multiplicities lower than or equal to 3. At high excita¬ 
tion energies the grand canonical ensemble is considered, 
as the number of partitions becomes too large. The choice 
of the partition is based on the statistical weights oc e^p 
in the first case, with the entropy Sp given by A. Botvina 
in [51] 



(18) 


and the average multiplicities of fragments of mass A and 
charge Z in the second case, 

{nA,z) = gA,zPA^/^elA(F.MT,v))-,A-.z], (19) 

Arp 

J-P 

These statistical weights and average multiplicities are 
based on several variable^ like the free energy Ea,z^ and 
a volume Vf related to the freeze-out volume V by V = 
Vb + V/, where Vq is the volume of the system with a reg¬ 
ular nuclear density. This volume Vf corresponds to the 
nucleus configuration where it is no longer influenced by 

^ TDHF: Time-Dependent Hartree-Fock 
^ QMD: Quantum Molecular Dynamics 
® SMM: Statistical Multifragmentation Model 
^ gA,z is the spin degeneracy factor, nA,z the number of 
fragments with mass A and charge Z in the partition, Aq (Eq) 
the mass (excitation energy) of the system, M the multiplicity 
of fragments, At^ the thermal nucleon wavelength, and y and 
n the chemical potentials. 
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the Coulomb forces and thus is no longer bound by the nu¬ 
clear forces. Then the fragments are dispersed by Coulomb 
repulsion from this freeze-out volume, after having been 
randomly positioned. 

Another approach has been studied in the de-excitation 
code Abla07 m- K.-H. Schmidt follows a practical way 
relied on experimental observations [65]. Even though, to 
determine the partition, it is not enough to take into ac¬ 
count the phase space only, knowledge of the dynamics 
are still poor. So, first of all, this multifragmentation pro¬ 
cess is not competing with the other de-excitation chan¬ 
nels and the remnant excited nucleus, derived from the 
intranuclear cascade, can undergo a multifragmentation 
if its temperature is above a threshold. In this case two 
options are offered: a fixed value for this temperature, or 
the use of a value depending on the mass of the excited 
nucleus [66]. The mass and charge of the emitted nuclei 
are obtained in two steps. The first step is the determina¬ 
tion of the amount of mass to be released (and thus the 
energy as well) to achieve a nucleus temperature below 
the multifragmentation threshold. This method is based 
on three points: i) When a unit mass is ejected, the nu¬ 
cleus loses energy ii) according to some experiments m 
this energy goes from 10 MeV for an excitation energy of 
2.9 A.MeV to 5 MeV for an energy of 11.8 A.MeV, and 
iii) the N/Z ratio of the initial nucleus is retained for the 
multifragmentation process, allowing to define the charge 
that must be emitted. The second step defines the masses 
Ai, charges Zi and excitation energies of the emitted nu¬ 
clei i. The masses are determined from a power law of the 
production cross-section observed experimentally 


da 


where the r parameter is taken in AblaOT linearly depen¬ 
dent of the excitation energy per nucleon. The masses are 
then defined iteratively until their sum matches exactly 
the ejected mass calculated in the first phase. Charges Zi 
related to the masses Ai are sampled in Gaussian distribu¬ 
tions in such a way that the mean charge Zmean respects 
the A/Z of the hot nucleus and that the width az de¬ 
pends on the temperature and on the symmetric term of 
the equation of state of the nuclear matter Cgym^ whose 
studies can be found in [68ll7Q] . 


a 


2 _ 
z — 


T 

Csym 


Finally, the excitation energy of the emitted nuclei {Z > 2) 
is based on the threshold temperature and on their den¬ 
sity parameter level, and their kinetic energy is due to 
the Coulomb repulsion and thermal agitation. For Z <2 
particles kinetic energy is drawn from a Maxwellian dis¬ 
tribution. 

As for the fission products, the various excited frag¬ 
ments from multifragmentation can in turn undergo evap¬ 
oration and/or fission for the heaviest ones and, for the 
lightest, in some models, a Fermi breakup. This latter pro¬ 
cess is described hereafter. 


2.2.5 Fermi breakup 

The Fermi breakup has been originally developed by E. 
Fermi to explain the multiple emission of tt in high ener- 
gie^ nuclear reactions m- He considered that the large 
number of possible states of the system leading to mul¬ 
tiple pion production justify a statistical approach. The 
probability of occurrence of a given channel was equal to 
its statistical weight (eq. [^, 


/ V \ 

statistical weight = ——^ pn{E)^ (20) 

\87r^h J 

with n the number of components (particles) in the final 
state, pn{E) the density of final states depending on the 
final total kinetic energy and V the volume of decay (often 
taken as the volume corresponding to normal nuclear den¬ 
sity). This method was then used for highly excited light 
nuclei [72], since the hypotheses were the same. Once a 
decay channel is selected, it remains to assign to each par¬ 
ticle its kinetic energy from the Coulomb repulsion, as in 
the case of multifragmentation. It is clear that the Fermi 
breakup and the statistical description of the multifrag¬ 
mentation are based on the same principle, a calculation 
of the probabilities of all possible partitions and selection 
by random sampling. Moreover, they involve both the vol¬ 
ume where the explosion happens, even if in the original 
Fermi case it was not a volume of freeze-out. Recently 
a study m, aiming at merging the Fermi-breakup model 
and the statistical multifragmentation model based on the 
microcanonical ensemble, showed that it was possible if 
the source volume was equivalent, which seems unfortu¬ 
nately not always the case. 


2.3 Transport codes 

Most people dealing with spallation reactions are not model 
developers, but users. For their projects involving such re¬ 
actions in massive targets, they need computational tools 
to get results. Moreover, they must also account for the 
particle transport in the given targets, as well as the slow¬ 
down, because the secondary particles can in turn induce 
other reactions. Thus the spallation codes must be imple¬ 
mented in a particle transport code, if the goal is not only 
to understand the elementary process, but also to offer 
a reliable tool to the scientific community. The transport 
codes, which manage the energy sector of spallation reac¬ 
tions must, to be complete, handle all the processes en¬ 
countered by a particle of energies ranging from a few GeV 
down to a few meV (for neutrons). For questions of knowl¬ 
edge of physics and for technical matters, two choices are 
possible to treat particle interactions: databases or mod¬ 
eling. Thus, regarding the interactions of low-energy par¬ 
ticles including neutrons, the choice of the database has 
been made. No model is able to reproduce the numer¬ 
ous resonances observed in the cross-sections, but, fortu¬ 
nately, the not too broad energy range (from a few MeV 

® at that time, a few hundreds MeV 
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to 20/150 MeV) allows experimental determination of the 
needed data, with the use of evaluations when necessary 
(see e.g. m)- In the spallation domain the situation is ex¬ 
actly the opposite. The energy range, from a few hundred 
MeV up to a few GeV, opens so many output channels 
that building a database fed by a battery of experiments 
is unimaginable, whereas our knowledge of the physics of 
spallation reactions provides an acceptable and now reli¬ 
able modeling of the interactions. This separation, or dif¬ 
ference, between low and high energies in the codes or 
systems of transport codes is illustrated in Fig. with the 
LAHET Code System (LCS) [75]. The high-energy part 
is handled by LAHET (Los Alamos High Energy Trans¬ 
port), based on the previous code HETC [TOhST] . itself an 
extension of NMTC [82| , whereas for the low-energy part 
the code HMCNP, modified version of MCNP [83] , is used. 

Another common point shared by the transport codes 
is the necessary reliability. They not only need to give the 
best results, which sometimes leads to a compromise be¬ 
tween sophisticated models and reasonable running time, 
but also cover an energy range broader and broader, in ad¬ 
dition to different operating systems evolving. The mainte¬ 
nance is crucial to keep a transport code alive. This aspect 
can have repercussions on the development of spallation 
models. 

As already said, a spallation model can be fully rec¬ 
ognized only if it can be used via a transport code, this 
means, if it is able to give results in all cases (projectile, en¬ 
ergy, target) and is reliable enough. Another point, rarely 
mentioned, is that its implementation in a transport code 
ensures in return the continued existence of such a model. 
Thus, for example, the INCL intranuclear cascade code 
has been implemented in several transport codes. Early 
in 2000 a collaboration with D. Prael allowed the use of 
INCL4.2 [27] and Abla [46] in the LAHET3.16 [76] code, 
and so in the LCS (see EigJ^. The goal was to be able 
to apply INCL combined with Abla on thick targets. This 
work done, and the first results being encouraging [84], the 
implementation in MCNPX2.5.0 (Monte Carlo N-Particle 
extended) was straighforward [85] . The particle transport 
codes handling spallation reaction have two origins: reac¬ 
tor physics and high-energy particle physics. Nowadays 
both families tend to join. Thus, willing to provide INCL 
to most of the people, it has been decide to implement it 
in GEANT4 (GEometry ANd Tracking) [86]. This code 
is a toolkit used by the high-energy physicists. The point 
was that only the codes written in C+-1- could be imple¬ 
mented. Then A. Boudard, co-developer of INCL, initiated 
the translation of the Eortran77 version of INCL in C+-1-. 
The work, started with P. Kaitaniemi m, has been con¬ 
tinued and finalized by D. Mancusi [88]. This huge work 
has been the opportunity to re-think the INCL code and 
has made its maintenance easier. In the mean time devel¬ 
opers of the PHITS (Particle and Heavy Ion Transport 
code System) transport code [89] decided to implement 
and use the INCL4.6 [28] version as their default intranu¬ 
clear cascade model. This example is mentioned to stress 
the point that the high quality of the physics model, if ab¬ 
solutely necessary, is unfortunately not enough to make a 



Eig. 8: The LAHET Code System (LCS) developed in Los 
Alamos, with the LAHET code for the high (spallation) 
part and the HMCNP code for the low energy. Eigure from 

m- 


computational code useful for the community. The coding 
part becomes then an essential activity. 

Coming back to the origin of the transport codes, the 
studies of meteorites clearly illustrate their convergence 
regarding the spallation reactions. Meteorites are inter¬ 
stellar bodies irradiated by the cosmic rays, mainly pro¬ 
tons, whose energy spectrum, going from a few MeV up 
to several TeV, peaks around 1 GeV, i.e. the spallation 
domain. We will come back later to such studies, but the 
energy range is so broad that people involved in this topic 
were interested in many types of codes. This is clearly seen 
by the use first of LCS by Reedy and Masarik |227j in 1994 
and then by the use of GEANT3[90] combined to MCNP 
m by Masarik and Beer m five years later. Einally 
Tab. gives an overview of the main particle transport 
codes dealing with spallation reactions with their charac¬ 
teristics. The ELUKA (ELUktuierende KAskade) [92] [93] 
and MARS (Math Analysis by Random Sampling) [941197] 
codes comes from the high energy physics. The former is 
developed mainly at CERN and the latter at Eermilab. 

Most of the time the lower energy limit is determined 
by the availability of data. Eor neutron projectile, databases 
exist for all nuclei below 20 MeV and sometimes also 
from 20 MeV to 150 MeV, whereas for light charged par¬ 
ticles the spallation model must be invoked down to 1 
MeV. Many more details on these transport codes can be 
found in the given references. These codes deal with en¬ 
ergy ranges broader than the spallation domain and so use 
either databases or other kinds of models. Some offer the 
possibility to mix data and codes when both exist (e.g. 
MCNPX). 
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Gode 

ING (Pre-equilibrium) 

De-excitation 

Projectiles 

Upper energy limil ® 

language 


Bertini (MPlvQl 


n, p 

TT 

3.5 GeV 

2.5 GeV 



Isabel (MPM) 

Dresneij^or Abla 

n, p 

TT 

d, t, ^He, a 

0.8 Ge\ 
1.0 Ge\ 
1.0 GeV/nu( 

c 

i 

zleorpl 

Eortran 90 

MGNPX2.7 

MGNP6 

INGL4.2 


n, p 

TT 

d, t, ^He, a 

-3 GeV “ 
-2.5 GeV 
—3 GeV/nucleon 



CEM03[IS2| + GEMflUH] 

n, p 

TT 

5 Ge\f 
2.5 Ge\ 



PHITS2. 64 

INGL4.6 

GEM 

n, p 

TT 

d, t, ^He, a 

3 Ge^l 

3 GeW 
3 GeV/nuc. 

leor0 

Eortran 77 


Bertini Intranuclear 

Gascade (+Preeq|^ 

Own Evaporatiorj^ (or GEM), 
Eission, 

Multifragmentatiorp] 
and Eermi breakup 
models 

n, p 

TT 

10 GeV 

10 GeV 


GEANT4 

Binary cascade (+Preeq.) 

n, p 

TT 

d, t, ^He, a 

10 GeV 

10 GeV 

few GeV/nucleon 

C++ 


INCL++ISSI 

or AblaV^ 

n, p 

TT 

d, t, ^He, a 

-3 GeV 
-3 GeV 

—3 GeV/nucleon 



PEANUT(GINC+preeq|^ 

Own Evaporation, 

n, p 

TT 

5 Ge\f 

5 GeV 

] - 


ELUKA 

r-QMD-2.^ 

Eission 

and Eermi breakup model^ 

d, t, ^He, a 

5 GeV/nucleorj^ 

Eortran 77 

MARS 

GEM03 

GEM 

n, p 

TT 

5 Ge\f 

5 GeV 

] - 

Eortran 77 


LAQGSM 

GEM 

d, t, '^He, a 

800 GeV/nucleon 



Table 1: Overview of the main particle transport codes dealing with spallation reactions. Sometimes, according 
to the codes, the multi choices (the physics models, their reliability in a given energy range according to a 
type of projectile) make it difficult to summarize in a table, thus for more details, readers are invited to look 
at the user’s manual of each code. While some codes allow the use of other projectiles than the nucleon, pion 
and composite particles up to the a particle, the list is limited to the projectiles that match the definition of 
spallation reaction given and addressed in this paper. 

® In some cases the limit can be changed by the user. 

^ The multistage multistep pre-equilibrium exciton model M- 

^ Above Bertini is used or FLUKA. It has to be mentioned that this is an old version of the FLUKA hadron generator, 
contained in FLUKA87[99]. That version has very little in common with the modern FLUKA [9^ [93] . 

^ Dresner is always combined with RAL|58] or ORNL[TOD] model for fission. 

® Like nucleons. 

^ Above the LAQGSM [lOlj model is invoked. 

^ Above the LAQGSM code is used. 

^ Above the JAM code uni] is used. 

^ Like nucleons. 

^ Above the JQMD code HSSl is used. 

^ Based on the work of Dostrovsky [106] . 

^ GEANT4 has its own pre-equilibrium model based on HOT]. 

^ Based on and [108] . 

Which is the Abla model previously mentioned, based on m and 

^ Above the Glauber-Gribov model insj is used. 

^ Based on the work of M. Blann [109] . 

^ This model handle nucleus-nucleus interaction below 5 GeV/nucleon [111] . 

Evaporation based on the Weisskopf-Ewing m formalism and fission on the Atchison model EH). 

^ Above the DPMJET model (based on [112] i is used. 

^ Above the LAQGSM code is used. 
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2.4 Conclusions on the modeling 

Most of the time the models describing the spallation reac¬ 
tions, or at least each phase, rely on the same basics, but 
the differences are in the details of the implementation, 
the way to do it and to put all the necessary mechanisms, 
ingredients, and only the necessary ones to be predictive 
and efficient. The knowledge of the spallation reactions is 
also driven by the need for simulations where they occur, 
in space or in accelerators, which gives another reason for 
getting the best description in terms of reliability and per¬ 
formance. The two next sections report on these aspects, 
with the validation of the spallation models, where their 
predictive powers prove their basic foundations, and their 
use and usefulness in various topics. 

3 Validations 

Model validation is an important step to know the value of 
a model. If everybody agrees with this preliminary remark, 
one must nevertheless, at first, define more precisely what 
is meant by validation. Validating can be understood in 
two slightly different ways depending on whether it is an 
end or a beginning. In other words, according to whether 
you are a model developer or a user, the goal of the vali¬ 
dation is not the same and nor is the way to do it. The de¬ 
veloper of a model sets goals to be achieved, he/she seeks 
to replicate such and such observable, and once he/she 
has built the model, he/she tests it on experimental data 
to judge its performance, to see if the mechanisms imple¬ 
mented fit the reality, and in this case the overall quali¬ 
tative approach is important. The user of a model wants 
to know before making calculations, which are the mod¬ 
els the most worthy of trust for the topic of interest and 
what are the margins of uncertainty. The approach is more 
targeted and quantitative. 

Trying to summarize the works on benchmarking of 
spallation models would be audacious, since all users did 
their own particular validations at a given time for a given 
project with one or several models, which makes the job 
too laborious and the information too specific. Therefore, 
in order to give an idea of the efforts regarding the vali¬ 
dation process of the spallation models and to tentatively 
draw some useful conclusions, in the following we draw 
up a brief history of the recent benchmarking campaigns, 
we address the work done by some transport codes teams, 
show the results of the most recent and most complete 
benchmark on elementary processes done under the aus¬ 
pices of the IAEA, and finally we mention some examples 
of validations in thick targets needing the use of transport 
codes. 

3.1 Overview of benchmark campaigns in the last two 
decades 

In the 1990s the Nuclear Energy Agency (NEA) of the 
OECD carried out a benchmark for three years, comparing 
codes in the field of spallation reactions. The goal was to 


get an inventory of modeling tools for spallation reactions 
in view of designing of accelerator-driven systems (ADS), 
dedicated to the transmutation of nuclear wastes. Two 
exercises were related to thin targets HH ITTK) . involv¬ 
ing microscopic spallation reactions, and one exercise was 
dedicated to thick targets thus including secondary 

reactions and therefore requiring the use of a transport 
code. The exercises on thin targets have attracted about 
twenty participants each. Chronologically the first HU] 
dealt with neutrons and protons emitted, while the sec¬ 
ond ESI dealt with the residual nuclei. A dozen models 
were implied in the thick target exercise. After these three 
years several conclusions could be drawn. The first obser¬ 
vation was that there was room for improvement in all 
models. Used nuclear masses and binding energies were 
not always the correct ones, a rough value being used 
instead of the measured value. To alleviate this prob¬ 
lem, details concerning level densities and shell effects had 
to be taken into account in order to correctly reproduce 
some reaction thresholds. Light composite particle yields 
were often badly reproduced and pre-equilibrium emission 
seemed to be needed to remove this deficiency. A more re¬ 
cent benchmark, discussed later in section [331 will show 
more precisely that this type of emission is necessary all 
along the intranuclear phase, but also that the role of the 
pre-equilibrium is not so clear. Modeling of light fragmen¬ 
tation products yields have also to be improved. The same 
is true for fission products as well. But in these cases the 
effort to remove deficiencies had to be done in the frame 
of de-excitation models. Einally models showed different 
behavior with low-energy projectiles. This is again due to 
the details in the models, as for the previous problems 
with reaction thresholds. A conclusion, not mentioned in 
the reports, implies that the models are globally good, 
and that it is time to turn to details. Moreover an evalu¬ 
ation of several models concerning a number of different 
observables with different targets and energies is always 
difficult. Eirst, very rarely, a model turns out to be bet¬ 
ter than another one for any of the observables; second, 
experimental data are not always sufficiently numerous to 
be representative of the whole set of possible reactions. 
Hopefully, the lack of experimental data will trigger new 
experiments that are expected to enlarge data sets. 

After these NEA benchmarks, knowledge and simula¬ 
tion of spallation reactions have shown significant progress 
in Europe, thanks to successive projects dealing with trans¬ 
mutation of nuclear waste and ADS. These projects in¬ 
cluded experiments intended to increase the amount of 
data and then to help improve the models of spallation. 
There first was the HINDAS |117j (High and Intermedi¬ 
ate Energy Nuclear Data for Accelerator-Driven System) 
project from 2000 to 2003, which benefited from a prelimi¬ 
nary previous work |118j . then NUDATRA |119j (NUclear 
DAta for TRAnsmutation) from 2005 to 2010 followed by 
ANDES |120j (Accurate Nuclear Data for Nuclear Energy 
Sustainability) from 2010 to 2013. 

In HINDAS, as in the other two projects, a report 
on the performance of the models was done at the be¬ 
ginning, firstly pointing out the difficulty in reproducing 
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some existing experimental data or the dispersion of pre¬ 
dictions of different models on observables where little or 
no experimental data exist, and secondly allowing defin¬ 
ing some specific domains to be improved. At the end 
of the project, models were tested especially on the new 
data obtained. If several models were tested at the begin¬ 
ning, like Bertini, Isabel (both implemented in the LA- 
HET Code System) and INCL for the intranuclear cascade 
and Dresner-Atchison and GEMINI for the de-excitation 
part, the spallation codes that were improved and com¬ 
bined were INCL for the intranuclear cascade, ending up 
in INCL4.2 [27], and Abla [46] for the de-excitation. This 
combination of codes has shown good progress in the sim¬ 
ulation of the spallation. It may be mentioned, among oth¬ 
ers, the ability of INCL4.2 to give the right reaction cross- 
section thanks to a diffuse nuclear surface, and a good 
reproduction of fission products for Abla through a new 
analytical approximation to the solution of the Eokker- 
Planck equation for the time-dependent fission width and 
a semi-phenomenological model m for the fission product 
yields. A coalescence model implemented in INCL, allow¬ 
ing cluster emission, showed already encouraging results in 
a previous version This idea will be applied later in 

an improved version of INCL |122j . This last point is also 
related to the excitation energy of the remnant, given by 
the intranuclear cascade model to the de-excitation model. 
This energy being a crucial input for the de-excitation, the 
intranuclear model must include the main channels in the 
right way to get it. It will be seen then that the pion emis¬ 
sion still have then to be improved in INCL. The new ex¬ 
perimental data obtained during this project will greatly 
help to validate the models, but also to show the way to 
improve them, and are still now considered as essential. 
We can cite the data on neutrons at SATURNE |123j and 
COSY |124j . the data on light charged particles at COSY 
|125j and the isotopic distributions for the residue produc¬ 
tion at GSI in inverse kinematics ISH HU niTi, comple¬ 
mentary to the excitation functions. More details can be 
found in references miiini. 

After HIND AS, INCL4.2 and Abla give generally rather 
good results, without free parameters, i.e. parameters de¬ 
fined once and for all, however both could still be im¬ 
proved. Thus some defects were to be corrected: produc¬ 
tion of light charged particles (up to a ), production of in¬ 
termediate mass fragments, called IME (with A < ^20), 
production of residues close to the target and handling of 
low-energy projectiles. This work was conducted within 
NUDATRA |119j . Although not mentioned in the list, the 
pion production in INCL has been improved by the works 
on nucleon and pion potentials [HHiiiisi, which is also a 
required step to get the right excitation energy. The estab¬ 
lishment of a process of coalesence allowed the production 
of composite particles up to alpha in a first step (INCL4.3 
|122j ) and up to nuclei with A<10 |13Qj later (but reduced 
in INCL4.5 [28] to A<8, for a question of computational 
time) and regarding Abla (become AblaOT) the evapora¬ 
tion of nucleons and a was extended to any type of par¬ 
ticles and nuclei in addition of a multifragmentation or 
breakup phase, open before the other de-excitation chan¬ 


nels and depending on a threshold temperature. Regard¬ 
ing the projectiles with energies below 200 MeV, great im¬ 
provements have been obtained with INCL by taking into 
account the Coulomb deflection and by handling correctly 
the soft collisions. This latter point is also mentioned as a 
good example demonstrating that a good prescription can 
prevent from improving a model if it is used outside of its 
original goal. In INCL the low-energy reactions, called soft 
collisions, are neglected, because they occurred often with 
no changes in the progress of the cascade, since most of 
the time they are Pauli-blocked. However this prescrip¬ 
tion must be removed, or at least modified, for the first 
collision of a low-energy projectile with a nucleon of the 
target nucleus, otherwise the number of events without 
interaction will be too high and the reaction cross-section 
significantly underestimated. Doing so, the INCL4.5 ver¬ 
sion was able then to take care correctly of the low-energy 
reactions. Einally, and unfortunately, the problem of the 
residues close to the target remained an open question. 
The experiments performed and used during this project 
were on light charged particle emission |131[ I132j and on 
the evolution of the fission cross-sections according to the 
projectile energy |133j . 

A little after the NUDATRA project it was decided 
to test the codes on a specific energy domain, i.e. in the 
range 150-600 MeV. This has been done in the ANDES 
project |120j . This energy range corresponds to the future 
prototype of the ADS Myrrha planned in Mol, Belgium. 
Data collected for the B ench mark of Spallation Models^ 
discussed later in section [33] were used with new experi¬ 
mental data for the energy in question. The details of this 
study are published in an ANDES report |134j . Eindings 
indicated good results, even if some defects were to be 
corrected and which tended to increase as the incident en¬ 
ergy decreased. Moreover, since on one hand the particles 
in a spallation target for ADS will be slowed down below 
150 MeV and the secondary low-energy particles can play 
a role, and, on the another hand, databases are missing 
for some types of particles in this energy range, the low 
energy interaction modeling will be also studied. Regard¬ 
ing intermediate mass fragment emission, if the INCL4.5- 
Abla07 gave rather good results, double differential spec¬ 
tra showed a high energy tail overestimated by INCL4.5. 
The coalescence model, where an outgoing particle ag¬ 
gregates with other nucleons, is based on a space-phase 
proximity criterion, and the use of a momentum cut-off 
seems to be needed as indicated by Eig. Improvement 
of clusters at low energy, together with residues close to 
the target, has been done considering the Eermi sea def¬ 
inition. In the cascade there are two types of nucleons: 
those which underwent a collision, the participants, and 
the others, the spectators. However, a participant that is 
slowed down below the Eermi energy is considered back 
to spectator. Actually the energy criterion is the Eermi 
energy plus a small quantity. In a first attempt this value 
was a constant value, 18 MeV, and the results were im¬ 
proved, but with a spurious hole in the neutron spectra. 
This was cured using a softening ”back-to-spectator” con¬ 
dition, i.e. requiring that the energy of the participant falls 
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Fig. 9: Double differential cross-sections for the produc¬ 
tion of ^Be for the reaction p(480 MeV)+Ag, as func¬ 
tion of the produced ^Be kinetic energy. Calculations 
with INCL4.6 (updated version of INCL4.5) combined to 
ABLA07 (red lines) are compared with the version includ¬ 
ing a momentum cut-off (blue lines) and with the experi¬ 
mental data (black symbols) |135j . Figure from |136j . 


below the emission threshold. This example showed once 
again that a careful implementation is required to take 
full advantage of a good idea. Regarding the secondary 
particles, INCL was able to emit clusters and handled in¬ 
teractions with particles up to a, but no special attention 
was paid to the low-energy light cluster collisions. Two 
aspects are important in this case: the scenarios (or en¬ 
trance channels) and some specific ingredients. The sce¬ 
narios that have to be taken into account are: formation 
of an (almost) undisturbed projectile fragment by non¬ 
interacting nucleons, formation of a compound nucleus 
(fusion) and cascade of interacting projectile nucleons; the 
two ingredients that become crucial at low energies are: 
Coulomb distortion and Q-values. This way the new ver¬ 
sion of INCL was made much more reliable in t his fie ld. 
More details are given in and in section |4.2.2l A 


last example of deficiency studied during ANDES was the 
overestimation of yields of some nuclei close to the tar¬ 
get and more precisely a typical case : the (p,2p) process. 
The overestimation of the cross-section could be improved 
by increasing the excitation energy, which facilitates the 
neutron evaporation. These (p,2p) events are probably pe¬ 
ripheral and the nucleons involved are the ones with the 
larger momenta, because usually nucleons are treated as 
classical particles, as far as their momentum is concerned, 
and sampled in a Fermi sphere distribution. Now, nucleons 
are quantum particles and the wave function of the less 
energetic nucleons in the nucleus can extend up to the 
surface. Thus, accounting for this quantum property, the 
removed peripheral nucleons could take away less energy, 
and so the excitation energy of the remnant nucleus could 


be increased. Encouraging results has been obtained in a 
recent study published by D. Mancusi |138j . where the 
quantum-mechanical features of the wave functions are 
mimicked by using a simple shell-model that allows deep 
orbitals to contribute to the surface energy content. 

In these three European projects the models, and more 
precisely INCL4 and Abla, were benchmarked and im¬ 
proved with some specific goals, but not used to design 
an experiment or a facility. In other words it was done 
mainly by the code developers. As mentioned in the intro¬ 
duction, one can also test models with the first aim not to 
improve them, but to use them in the best way, that is to 
choose the most suitable model, knowing its uncertainties 
concerning a studied subject. Eor example we can cite the 
work done within the EURISOL-DS (EURopean Isotope 
Separation On-Line Radioactive Ion Beam Eacility - De¬ 
sign Study) project, from 2005 to 2009. The project aimed 
at designing a facility delivering radioactive ion beams us¬ 
ing two types of spallation targets: one producing the ions 
directly from spallation reactions and the other using the 
spallation neutrons to induce fission in a uranium target. 
Several types of materials were to be studied and at first 
a benchmark of the tools, spallation codes, was done. As 
already said, there is no use here to discuss all results, but 
it is worth saying that for this project three of the five 
transport codes listed in Tab. were studied as well as all 
model combinations in MCNPX2.5.0. The reader will find 
information, concerning MCNPX2.5.0, in the references 
[numi, where it was concluded that the two model 
combinations Isabel-Abla and INCL-Abla, as well as the 
CEM2k model, were the three best options for this pur¬ 
pose. Regarding ELUKA, it was first benchmarked on the 
elementary residue production cross-sections, reproducing 
the experimental data within a factor of two in most of the 
cases |142j . like the best options of MCNPX2.5.0. Then it 
was compared to neutron multiplicities |143j with more or 
less good results and finally validated on neutron attenu¬ 
ation in shielding materials nm where the transport of 
the particles is important. In |145j a comparison of PHITS 
with MCNPX2.5.(f]for some ion-ion interactions was per¬ 
formed. The question of considering these latter reactions 
as spallation reactions is still open. The numerous possibil¬ 
ities concerning the target materials forced researchers to 
search for a large set of experimental data, which will be 
efficiently reused for the Benchmark of Spallation Mod¬ 
els. However, the benchmark of elementary processes in 
the transport codes within EURISOL-DS was obviously 
focussed on the projectile energy, i.e. around 1 GeV, and 
more benchmarks can be found on the websites of the de¬ 
velopers of these transport codes as discussed in the next 
section. 

Regarding the experimental data, used for such vali¬ 
dations, great efforts have been made to collect and offer 
high quality data. We can cite for the elementary reactions 
the Experimental Nuclear Reaction Data base EXEOR 
|146j . This database originally dedicated to low-energy 
neutron cross-sections has been extended to light charged 
particles and to higher energies, say around 1 GeV, but 

^ Using the LAQGSM model in this case 
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sometimes beyond. EXFOR aims at providing experimen¬ 
tal data with all necessary information, when available, to 
make them usable. EXFOR is a network of Nuclear Re¬ 
action Data Centres (NRDCp^ coordinated by the Inter¬ 
national Atomic Energy Agency (IAEA). Exchange of in¬ 
formation between developers, users and experimentalists 
are promoted by the Working Party on international nu¬ 
clear data Evaluation Co-operation (WPEC) established 
by the NEA Data Bank in close co-operation with the Nu¬ 
clear Data Section of the IAEA. For example a subgroup 
dealing with Quality Improvement of the EXFOR database 
gave the opportunity to improve the data related to 
spallation reactions (some errors were corrected) and the 
way to find thenj^ On the side of thick, massive or com¬ 
plex targets, each project is a new one, however some typ¬ 
ical cases exist and a joint project was launched in 1996 
by the NEA Data Bank and ORNL/RSICC (Oak Ridge 
National Laboratory, Radiation Safety Information Com¬ 
putational Center): SINBAD (Shielding Integral Bench¬ 
mark Archive Database). The objective was to produce a 
database containing an internationally established set of 
radiation shielding and dosimetry data. For example, the 
experimental data used within EURISOL-DS for the val¬ 
idation of FLUKA on neutron attenuation came from the 
HIM AC experiments included in SINBAD. More informa¬ 
tion can be found in 

Validation is a long and complex process, starting with 
numerous benchmarks in laboratories done privately or 
within collaborations on international projects followed by 
the work of the transport code developers, and, from time 
to time, with international benchmarks as the ones of NEA 
and IAEA making the point on the spallation reaction 
modeling. The following sections give more details on some 
examples to show the difficulties and the outcomes of the 
different types of benchmarks. 


3.2 Benchmarks within transport codes 

Particle transport codes dealing with spallation reactions 
are tools supposed to work with all kinds of material (or 
target), all types of projectiles, with a module handling 
the lower energy domain when the particles are slowed 
down to a threshold depending on the spallation models. 
Moreover most of these codes also manage higher energy 
reactions. Thus a benchmarking or validation activity can¬ 
not aim at covering all possible studies, as it is written in 
the MCNPX2.6.0 User’s Manual |148j : All eomputer sim¬ 
ulation eodes must he validated for speeifie uses, and the 
needs of one projeet may not overlap eompletely with the 
needs of other projeets. It is the responsibility of the user 
to ensure that his or her needs are adequately identified, 
and that henehmarking aetivities are performed to aseer- 
tain how aeeurately the eode will perform. The henehmark¬ 
ing done for eode developments for the MCNPX sponsors 

https://www-nds.iaea.org/nrdc/ 

Some data come from experiments done in inverse kine¬ 
matics, this means that the heavy nucleus is the beam and the 
light particle the target. These data must be also found when 
considering the light particle as the beam. 


may or may not be adequate for the needs of the user^s 
partieular program. We make our henehmarking efforts 
publie as they are eompleted, but the user must also de¬ 
velop a rigorous henehmarking program for his/her own 
applieation. Nevertheless great efforts have been carried 
out by the transport code developers to give indications 
on the reliability of their code, often motivated by par¬ 
ticular applications or current topics. We give hereafter 
some examples of the work performed by three of them: 
PHITS, MCNPX/MCNP6 and GEANT4. 


Around forty references related to benchmark studies 
can be found on the website of the PHITS code ms]. 
PHITS covers a large energy scale from say 0 to 200 GeV 
with not only light hadrons as projectiles, but also heavy 
ions, leptons and gammas, and the work is almost always 
motivated by a specific goal, e.g. safety, medicine, cos¬ 
mic rays, accelerators. Thus, if more than ten papers deal 
with benchmarks with thick target (shieldings studies) 
and more or less the same amount for heavy ion reactions, 
the number of papers on the validation of the elementary 
processes of the spallation models in PHITS is far from 
being the majority. However these studies deal with diffi¬ 
cult points and experiments are sometimes carried out to 
provide the needed experimental data, which are then use¬ 
ful for the whole community. Hereafter some examples are 
cited, first to mention what has been done by the PHITS 
developers, but also to address some points raised by these 
benchmarks and discussed later. It must be added that the 
PHITS team was involved in the Benehmark of Spallation 


Models developed in section 3.3 


Several spallation models are available in the PHITS 
code and benchmarks are often used to know which one 
is the best in a particular case. It is not the place here 
to summarize the history which ends up in the choice of 
the default options or default models in PHITS, and all 
information is in However three interesting points, 

at least, can be drawn from this validation work. The first 
one concerns the libraries, the second one the neutron- 
induced spallation reactions and the last one the pions 
and resonances. Below 20 MeV or 150 MeV it is common 
to use data libraries instead of models, the former ones 
being much faster and validated. Beyond, the problem is 
the large energy range and the numerous output channels, 
which makes the building of such libraries tedious. Never¬ 
theless this topic is developed later in section [4.1.1[ and we 
can mention here that, at least in two papers, PHITS de¬ 
velopers validated their own JENDL-HE libraries. Briefly 
the GNASH |150j code has been used below 150 MeV 
and the JQMD model above, and up to 3 GeV, to con¬ 
struct this library. In the paper dealing with forward neu¬ 
tron emission from the reaction p(100-200 MeV)+G |151j 
the JENDL-HE library shows its ability to reproduce the 
Quasi-Free peak contrary to the spallation models. This 
is due to the choice of the right and dedicated models to 
build the library in the specific energy ranges. However 
this asset can become a drawback: a gap at the frontiers 
of the models used. This appears in the neutron activa¬ 
tion cross-sections study where a discontinuity is 

seen around 250 MeV. The experimental data used for this 
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study were measured on purpose, because, while below 20 
MeV numerous experimental data are available, they are 
scarce between 20 MeV and 200 MeV and nothing exists 
above. Thus several experiments had been carried out at 
the Research Center for Nuclear Physics (RCNP) at Os¬ 
aka using a quasi-monoenergetic neutron beam of around 
270 MeV and around 380 MeV (references in |152j h This 
kind of data is strongly needed in particular for extrater¬ 
restrial matter irradiation by cosmic rays, and some other 
techniques have been developed to extract such informa¬ 
tion |153j . The spallation model used by PHITS in this 
study |152j is INCL4.6 [28] and the results are reasonable, 
with obviously no discontinuity in energy. This specific 
point on excitation functions with neutron projectiles is 
discussed in section |4.1.2| A third point, discussed in an¬ 
other validation of PHITS on neutron production via the 
pion interactions |154j . is on the better way to take into 
account the resonances. The place where the Z\(1232) is 
produced in the nucleus is invoked to explain the produc¬ 
tion of negative and positive pions in iron and lead, as well 
as the nucleon resonances taken into account in JQMD to 
understand the differences between the two types of pions. 
This last point is an example of what must be included 
to get the best results in a simpler and faster way. The 
production of pions around GeV, like in this paper, and 
at higher energy can be handled in two ways: either by the 
decay of the resonances considered as participants of the 
intranuclear cascade or forgetting this intermediate step 
and taking into account the direct production of pions. 
So, for example, the first way is implemented in JQMD, 
whereas the INCL code used the second method. This way 
to proceed is explained in more details in sections [5.1 [and 


Like PHITS, MCNPX/MCNP6 put in its website the 
references of the works done about validations. Around 
fifty papers are cited and the half dedicated to the spal¬ 
lation reactions. The rest is mainly on low energy neu¬ 
tron and criticality. The validations of the spallation mod¬ 
els in MCNPX/MCNP6 have been carried out most of 
the time by S. Mashnik, focussing on the GEM code, 
but also LAQGSM at higher energies, with comparison to 
experimental data and other codes implemented in MG- 
NPX/MGNP6 or not. So, through these benchmarks, we 
can follow the development of the GEM code, from GEM95 
|155j to GEM03.03 [29| via GEM2k |156j . There is no 
doubt that this continuous work of validation to test and 
to know the points to improve is the reason why in MCNP6 
the CEM03.03 has become the default option, replacing 
the famous Bertini-Dresner combination. Among all these 
references HSU, we can cite the one on validation and 
verification of MCNP6 |158j , where eighteen tests are pre¬ 
sented in details with the input files to help the users. 
Light and heavy targets with mainly protons, but also 
neutrons and pions as projectiles are considered in a wide 
energy range. Three tests are dedicated to gamma-induced 
reactions around GeV. Photonuclear reactions are not usu¬ 
ally considered as spallation reactions, but this interest¬ 
ing topic is addressed in the last section [^ The fifteen 
other benchmarks were mainly dedicated to the elemen¬ 


tary processes with thirteen tests, two tests dealing with 
thick targets. In addition, considering the spallation reac¬ 
tions as the interaction of a light particle with a nucleus, 
we can mention two other tests described in another report 
dedicated to the LAQGSM code One on backward 

neutron emission from a thick lead target, a difficult task 
for the spallation models, and another one induced by al¬ 
phas, not handled by GEM03.03. Doing so the codes have 
been verified, i.e. can work with no bugs, and validated, 
i.e. give reasonable results and can be considered reliable. 

The website of GEANT4 offers also results of bench¬ 
marking works. However, some differences exist with the 
two previous transport codes. If in the three cases the 
references of the experimental data used are given, only 
comments on the results are available in the GEANT4 case 
(one article can nevertheless be mentioned [160]), prob¬ 
ably because, in return, an impressive number of plots 
from the calculation results have been obtained and are 
put in the website, which allows anybody to make its own 
analysis. We see clearly with the GEANT4 validation pro¬ 
cedure, compared to the other two, that a compromise 
has to be made between the number of tests and the de¬ 
tail of the analysis. Once again the validations can be 
undertaken with different goals: user-oriented, developer- 
oriented or in-between. GEANT4 seems the most user- 
oriented regarding the amount of information provided, 
whereas PHITS and MCNPX/MCNP6 give maybe less 
resultbut with most of the time a critical analysis. 
There is however one difficulty that each has to face, which 
is the regular updating of such benchmarks. Nevertheless 
the work done and the information made available by the 
transport codes concerning validation is a great database 
for the user as well as for the developer. 

Einally we can mention two papers dealing with an 
inter-comparison of transport codes, in particular those 
in Tab. but at the frontiers of the spallation reactions. 
One benchmarked heavy ion reactions and the other 
one high energy reactions |162j . i.e. the lowest projectile 
energy being 12 GeV/c. In those works the number of 
reactions were limited to draw clear conclusions. 


3.3 IAEA Benchmark of Spallation Models 

Around fifteen years after the benchmarks carried out by 
the NEA and followed by great improvements in the spal¬ 
lation reactions modeling, three institutes, CEA (Erance), 
Soreq (Israel) and EZJ (Germany), decided to launch an¬ 
other validation exercise under the auspices of the IAEA. 
In order to get the most discerning conclusions the goal 
was to focus on the elementary processes only. Neverthe¬ 
less, trying to make an assessment of the spallation models 
used around the world can be an ambitious task for many 
reasons: 

• scanning all areas of spallation reactions (type of pro¬ 
jectile, projectile energy, type of nucleus target, ob¬ 
servable) 

difficult to compare, because the cases studied are not al¬ 
ways the same ones 
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• involving the most of existing models, 

• asking all participants to perform all calculations with 
identical output formats, 

• creating the tools for the analysis of results, 

• making the most comprehensive analysis, 

• making this exercise valuable to developers and users, 
and 

• making available the information and tools used for 
future next exercises. 

This IAEA Benchmark of Spallation Models is developed 
hereafter for several reasons. First, up to now, it is the 
most complete benchmark in this field, and second, if some 
papers have been written to summarize the results ob¬ 
tained some information provided during the 

exercise can still be analyzed and some points are un¬ 
derlined here. The author apologizes in advance for the 
resolution of the plots shown in this section, most of the 
time taken from the website of the benchmark. The goal 
of these plots is here often to show the trend, not to dis¬ 
cuss a particular model, then if the reader requires more 
information, all information is available on the website. 

A first workshop was held in Trieste in February 2008 
to define the content of this benchmark. It brought to¬ 
gether developers presenting their models and experimen¬ 
talists who made a review of existing data. The first result 
of this workshop is a report m with, in particular, inter¬ 
esting descriptions and details of many models. Another 
important point was the compromise on the experimen¬ 
tal data that should be used for the comparisons. The 
aim was to cover the broadest spectrum of spallation re¬ 
actions taking into account two constraints: the availabil¬ 
ity of data and the amount of work that would lead to 
the participants and the organizers. Taking into account 
these objectives and constraints, it was decided to con¬ 
sider as projectiles only nucleons, mainly protons with in¬ 
cident energy from 20 MeV to 3000 MeV, but, as targets, 
heavy, intermediate and light nuclei. The observables were 
spectra of particles (neutrons, protons, deuterons, tritons, 
^He, a and pions), neutron multiplicities, isotopic (mass 
and charge) distributions of the residual nuclei, and iso¬ 
tope production function of the energy of the projectile 
(excitation functions). 

Table lists the experimental data sets. These data 
were then also put in the ASCII format to be easy to 
use during this benchmark, but also by anyone, since the 
files have been made available via a dedicated website 
(http://www-nds.iaea.org/spallations/). Each of these files 
has been associated with another one summarizing the 
necessary information on the experimental conditions. As 
regards the energy range under consideration, the domain 
below say 150 MeV has been validated also, because, when 
used in a transport code, the models are called when no 
database exists, and especially at low energies. Also, it is 
interesting to know at least the behavior of the models in 
these domains, where they are not supposed to give good 
results. 

Any person involved had to commit to do the calcula¬ 
tions for all the reactions mentioned in the exercise, in or¬ 
der to get a complete and consistent computational results 


database. Regarding the excitation functions it was agreed 
that the participants gave the direct (or independent) pro¬ 
duction cross-sections, while the experimental data also 
included cumulative productions, which take into account 
the radioactive progenitors. Actually only one model for 
the calculation of the cumulative cross-sections was used 
by the organizers to avoid bias in the analysis of these ob¬ 
servables. Details on this cumulative model can be found 
in |194j . Table 1^ shows the models or model combinations 
involved in this benchmark and the calculation results 
provided (in green) or not (in red). The first result was 
twofold: the most widely used spallation models (except 
the one in the FLUKA transport code) were involved and 
a complete picture of the spallation reactions modeling 
could be drawn. 

In order to analyze the calculation results a huge num¬ 
ber of figures (several thousands) have been built. Most of 
the figures are the results of each model compared to the 
experimental data, but the results of all models were also 
plotted on a common figure, for a given experimental data 
set, to make a global inter-comparison. We will see here¬ 
after that this type of figure allows the trends to be seen, 
if the experimental data are collectively well reproduced 
or not, if the models give more or less the same results, 
with the same quality or shortcomings, or if the results 
are very different. In addition to the usual figures. Figures 
of Merit (FoM), or deviation factors, were also drawn for 
each model and observable. Those tools deserve some at¬ 
tention, because they are symptomatic of the goal of the 
validation. Classical figures comparing the calculation re¬ 
sults with experimental data of an observable can be seen 
as a global and qualitative approach. The FoM are quan¬ 
titative, but focused on one characteristic. The choice of 
the tools is driven by the aim of the validation. The FoM 
that were used are the following factors: 



where is the measured (calculated) cross-section, 

the experimental uncertainty and N the number of 
experimental data points. Each of these coefficients shows 
a particular characteristic of the calculation result. R is 
the crude mean ratio, and is not symmetrical, since it is 
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Projectile 

Target 

Energy (MeV) 

Observable 

Reference 

n 


65 

DDXS - neutron 

m 

P 

2D5p^^ 

63 

DDXS - neutron 

nMi 

P 

natpb 

256 

DDXS - neutron 

mii 

P 

natpe 

800 

DDXS - neutron 

um 

P 

natpg natpp 

800, 1200, 1600 

DDXS - neutron 

Multiplicity - neutron 

mi 

P 

^®^Pe ’^“^Pb 

3000 

DDXS - neutron 

UM 

P 


1200 

Multiplicity distribution - neutron 

um 

P 

noTpp 

1200 

Multiplicity distribution - neutron 

HZQ] 

n 

2U9Bi 

542 

DDXS - proton, deuton, triton 

mi 

P 

'’“Fe, 

62 

DDXS - proton, deuton, triton, ^He, a 

mi 

P 

2nsp^^ 

63 

DDXS - proton, deuton, triton, ^He, a 

nssi 

P 

^'Al, '“"Au 

160 

DDXS- a 

m 

P 

nat^i 

175 

DDXS - proton, deuton, triton, ^He, a 

um 

P 


175 

DDXS - proton 

11151 

P 

2D5p^^ 

800 

DDXS - proton 

nmnzzi 

P 


1200 

DDXS - proton, deuton, triton, ^He, a 

um 

P 

^^^Au 

1200 

DDXS - proton, deuton, triton, ^He, a 

m 

P 

^^^Au 

2500 

DDXS - proton, deuton, triton, ^He, a 

|121||131| 

P 

"“*C, ^‘'Al, ’"“*Cu and "“‘Pb 

730 

DDXS - 7r± 

nza 

P 

^’Ai 

2205 

DDXS - TT- 

mi 

P 


300 

Isotopic distributior ® 

mi 

P 


1000 

Isotopic distribution 

|1811I182| 

P 

2D5p^^ 

500 

Isotopic distribution 

um 

P 

2nsp^^ 

1000 

Isotopic distribution 

m 

P 

238 p 

1000 

Isotopic distribution 

imiinziiiHaiiHS] 

P 



Excitation function 

111511186H190I 

P 

natpb 


Excitation function 

|191H19;i| 


Table 2: Experimental data sets considered in the IAEA Benchmark of Spallation Models. 


Experiments done in inverse kinematics 


Code 

Representative 

Neutron 

DDXS 

Neutron 

Multiplicity 

p, d, t, 
®He, a 

Pions 

Residue 

distributions 

Excitation 

functions 

cem03-02 

A. Sierk, S. Mashnik, 

K. Gudima, M. Baznat 







cem03-03 

K. Gudima, M. Baznat, 
A. Sierk, S. Mashnik 







cascade04 

H. Kumawat 







cascade-asf 

A. Konobeyev 







cascadex 

Yu. Korovin 







phits-bertini 

phits-jam 

phits-jqmd 

N. Matsuda 







Geant4-bertini 

Geant4-binary 

Geant4 Hadronic 
Group 







isabel-abla07 

isabel-gemini-h+ 

isabel-smm 

D. Mancusi 







incl4.5-abla07 

incl4.5-gemini-|-+ 

incl4.5-smm 

D. Mancusi 







MCNPX 

(Bertini-Dresner) 

F. Gallmeier 


No multiplicity 
distributions 






Table 3: Spallation models or model combinations involved in the IAEA Benchmark of Spallation Models. More 
details on the models (versions, ingredients, ...) at https://www-nds.iaea.org/spallations/2010ws/ 
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Acceptance band [eval/x ; eval*x] 

Points 

x=5 

1 

x=3 

2 

x=2 

3 

x=1.4 

4 


Quality 

Points 

Good 

2 

Moderately good, minor problems 

1 

Moderately bad, particular problems 

-1 


Unacceptably bad, systematically wrong -2 


Table 4: Ratings used to analyze the benchmark results: 
upper part for light charged particles and lower part for 
neutrons and residues. 


dominated by the overestimation, H is the mean F is 
a symmetric version of R, i.e. overestimation and underes¬ 
timate treated on an equal footing, S is similar to F, but 
taking into account the weight of the experimental errors, 
M, called intrinsic discrepancy Una, can be seen as the 
F factor, but assessing the shape of the distribution, and 
finally F^,, which allows to get a rough idea of the accuracy 
of a model for a given calculation. 

Two types of analysis have been carried out and pre¬ 
sented during the ” Second Advanced Workshop on Model 
Codes for Spallation Reactions” at Saclay (France) in 2010. 
On one hand each participant was asked to explain the 
successes and deficiencies of their models and, on the other 
hand, a global analysis showed the main trends and a sta¬ 
tus of the spallation modeling. This last task has been 
divided in three subtasks, neutron, light charged particles 
and residues, taken care of respectively by J.-C. David 
(CEA-France), F. Gallmeier (ORNL-USA) and R. Michel 
(ZSR-Germany). A summary of this work can be read in 
[muM], but the great deal of information built up in 
the website dedicated to this benchmark can be still used 
to study more details. Here is a picture we can paint of 
the spallation modeling based on this benchmark. 

Regarding the global analysis, two methods were used 
to rate the models. The first one was suggested and used 
by F. Gallmeier for the light charged particles. It is based 
on a quantitative procedure. The more the calculation re¬ 
sults deviate from the experimental data, the lower is its 
rating. The second one was suggested by R. Michel and 
based on a more qualitative method. The ratings for both 
methods are given in Table Obviously the plots, and 
sometimes the FoM, were studied to go into the details 
(ingredients, mechanisms, ...). 

The data studied for the neutrons are mainly the dou¬ 
ble differential cross-sections, to which the rating proce¬ 
dure was applied, but also the multiplicities and multiplic¬ 
ity distributions. In order to separate the different mecha¬ 
nisms, each spectrum was divided into energy (0-20 MeV 
related to the de-excitation step; 20-150 MeV for the inter¬ 
mediate phase sometimes called pre-equilibrium; beyond 


150 MeV for the intranuclear cascade) and into angle (for¬ 
ward angles (< 45°); intermediate angles (45° < 0 < 135°), 
and backward angles (> 135°)). Figure 10 shows the re¬ 
sult of ratings averaged over these regions and over all 
the reactions studied. One can clearly see that the disper¬ 
sion between the models is small and the overall results 
are correct. However, a more detailed study of the spectra 
can point to deficiencies and differences. 



Mean per region 


1 

■ 




1 

1 

■ 

■ 

■ 

1 

1 

1 

1 

1 

1 

1 

1 

1 
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Fig. 10: Rating results obtained using the method given 
in the lower part of Table for neutron double differential 
cross-sections. 


In the low-energy part the models give good results. 
This is because evaporation of neutrons is the preferred 
route of the de-excitation and the formalism of evapo¬ 
ration is well known. Nevertheless two features may be 
noted. On the one hand the border between the neutrons 
from the intranuclear cascade and the neutrons from the 
de-excitation is very difficult to locate precisely, thus to¬ 
wards 20 MeV the differences between the models are more 
pronounced, and on the other hand, the other side of the 
spectrum to very low energies, the models can differ a 
lot. This last point is possibly linked to the capture cross- 
section of neutrons used in evaporation models. The val¬ 
ues used for these cross-sections rest most of the time on 
experimental measurements obviously made with cold nu¬ 
clei, while here the nuclei are hot, therefore the energy 
levels are different and the probability of low-energy neu¬ 
tron capture modified. 

The portion of the spectrum between 20 MeV and 150 
MeV is linked to the end of the intranuclear cascade. Some 
spallation models add an intermediate step between the 
INC and the de-excitation, the pre-equilibrium, believing 
that the remnant nucleus from the cascade does not have 
all the characteristics of a compound nucleus. Figure pT] 
shows clearly below 100 MeV two groups of model. If we 
exclude the two codes JQMD (from the PHITS transport 
code) and CASCADEX, these two groups are character¬ 
ized by having or not a pre-equilibrium phase. One way to 
highlight the usefulness or not of a pre-equilibrium model 
is to look at, for example, the neutron spectrum in the re¬ 
action p(63MeV) + Pb, where the incident energy is such 
that models with pre-equilibrium should give better re¬ 
sults. On the upper part of Eigure[^are shown the results 
of the models with a pre-equilibrium step and on the lower 
part the models without pre-equilibrium. While the re- 
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Neutron energy (MeV) 


Fig. 11: Double differential cross-section production of 
neutrons at an angle of 85° from the reaction p (IbOOMeV) 
+ Fe. This plot focuses on the intermediate energy part 
^0 MeV - 150 MeV). The seventeen models (listed in Tab. 
1^ are plotted and the experimental data taken from |123| . 


suits, for three angles, of the models with pre-equilibrium 
are less scattered and give best results for the backward 
angle, the models without pre-equilibrium are generally 
better at the forward angle. The situation is more com¬ 
plicated for the transverse direction. The diversity of the 
results for models without pre-equilibrium is due to the 
special treatment of this low incident energy, outside the 
field of the spallation reactions. Besides, a model with¬ 
out pre-equilibrium, INCL4.5, gives interesting results if 
we ignore the spurious fall at ^10 MeV. This defect in 
INCL4.5 was corrected in the following version INCL4.6 
leading to good results, as shown in Fig{^. From all this 
it appears that the usefulness of a pre-equilibrium model 
is not obvious. 


The high-energy part of the neutron spectrum is very 
anisotropic, thus one must distinguish the forward angles 
(< 45°- Fig. 14), the intermediate angles (45 < 0 < 135°- 
Fig. 15) and the backward angles (> 135°- Fig. f^. 


The forward angles are characterized by two peaks, 
which disappear gradually as the angle increases (Fig.p!4|). 
The quasi-elastic peak is located at energies close to that 
of the incident projectile, it is the consequence of the direct 
interaction of the projectile with a neutron in the nucleus 
and thus is sensitive to the nucleon-nucleon cross-sections 
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Fig. 12: Double differential cross-section production of 
neutrons at three angles (24° on the upper part, 80° in 
the middle and 120° on the lower part) from the reaction 
p(63MeV) + ^^^Pb. The left part deals with the models 
including a pre-equilibrium step and the right part the 
models without pre-equilibrium. Experimental data taken 
from |166j . 


used, as well as to the spatial distribution of the nucleons 
in the nucleus. While the peak position is often correct, its 
height is a problem for many models. The quasi-inelastic 
peak is located slightly below the quasi-elastic peak and 
is significantly broader. This peak is due to the decay of 
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TfMeV) 

Fig. 13: Comparison of INCL4.5 (blue) and INCL4.6 
(red) predictions for the neutron double differential cross- 
sections in p(63MeV) + Pb collisions. For the sake of 
clarity, the various spectra have been multiplied by 10^, 
10^, 10^, etc., for increasing angles, starting at the lowest 
angle. Experimental data taken from |166j . Figure from 
[28] , Copyright (2014) by the American Physical Society. 


the Z\(1232) resonance and does not seem always easy to 
reproduce at all angles. At 25° some models still exhibit 
a peak in contrary to the experimental data and at 0° 
some differences of one order of magnitude can exist on 
its height between models and measurements. However 
it is worthwhile to note that at 0° the simulations are 
difficult, because a compromise has to be found between 
the optimum aperture angle and the computational time. 

The intermediate angles are the part of the high en¬ 
ergy spectrum reproduced the best by all models. Angu¬ 
lar distributions of the FoM M (formula page 20) are 
shown in Fig. 15 for the reaction p(800 MeV) -j- Fe and 
one clearly sees that the models reproduce much better 
the transverse angles (M is close to 0) than the forward 
and backward angles. 

One example of spectra of neutrons emitted at high 
energies in the backward direction is given in Fig. If 
the shape of the distribution is fairly well reproduced, the 
models differ by the values of the cross-sections first, and, 
secondly, towards the ’’lower” energies some models devi¬ 


ate from the experimental data. This last point is probably 
due to the multiple scattering. Other data, not shown here, 
point out that these difficulties increase with the angle. 

The neutron spectra can be used to address the high 
energy limit of the spallation models. This is illustrated by 
Fig. where neutron spectra at the same angles (close 
to 90^ are plotted with the same target, but with two 
incident energies: 1.2 GeV and 3.0 GeV. This kinematics 
region of neutron production is, as previously stated, well 
reproduced by the models (see the left part), but when the 
incident energy increases and reaches the value of 3.0 GeV, 
the calculation results are scattered and less good. The 
reliability of the spallation models above 3.0 GeV depends 
strongly on the models. The point, to extend the models to 
higher energies, is to take carefully into account the needed 
channels in the intranuclear cascade (resonances?) and in 
the de-excitation phase (multifragmentation?). This point 
is discussed in section [5T1 

To conclude on the neutron production Figure [Ts] shows 
six neutron multiplicities, three incident energies (800 MeV, 
1.2 GeV and 1.6 GeV) for low-energy neutrons with an 
iron target and the same three incident energies for high- 
energy neutrons with a lead target. The same conclusions 
as those from the double differential spectra can be drawn. 
So, if in the case of iron some rare calculation results are 
within a factor of 2, whereas many models are within 30% 
or better, in the case of lead nearly all points are in the 
error bars or very near. Actually the differences must prob¬ 
ably come from the neutron energies and not from the tar¬ 
get. The low-energy neutrons come from the de-excitation 
phase and are in competition with all other de-excitation 
channels. The multiplicities are most of the time overesti¬ 
mated, which could indicate the lack of other channels. 

The light charged particles considered during this bench¬ 
mark were the proton, the deuteron, the triton, the ^He 
and the a. The rating procedure, explained in the up¬ 
per part of Tab.[^ was applied to each data set (no global 
rating). Fig. [l^ shows two examples, which allow the main 
conclusions to be drawn on the modeling of the light charged 
particle emission. While some models seem to be generally 
more efficient than others, the results depend very much 
on the reaction and on the given emitted particles. Thus, 
for the reaction n(542 MeV) + Bi (upper part of Fig. 


19) protons are reasonably reproduced by many models. 


which is not the case at higher energies, as seen with the 
reaction p(2500 MeV) -f Au (lower part of Fig. 19). Al¬ 
though the emitted particle energy range is not exactly 
the same in both reactions, data at 542 MeV not cover¬ 
ing the low-energy (evaporation) region, it is clear that 
the incident energy plays a role. The case of the compos¬ 
ite particles (d, t, ^He and a) clearly divides the models. 
There are those that emit the charged particles through¬ 
out the spallation reaction and those which do not have 
a suitable mechanism for high energy emission (intranu¬ 
clear cascade). Every developer agreed on the need for a 
coalescence model in the ING to be able to emit cluster 
particles. 
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Fig. 14: Double differential cross-section production of neutrons at three forward angles (0° on the left, 
10° in the middle and 25° on the right) from the reaction p (1200MeV) + Fe. This plot focuses on the high 
energy part (above 100 MeV). The seventeen models (listed in Tab. ^ are plotted and the experimental 
data taken from |123j . 
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Fig. 15: Angular distributions of the FoM M related to the high-energy (>150 MeV) neutrons produced 
in the reaction p (SOOMeV) + Fe. The seventeen models (listed in Tab.l^ are plotted and the experimental 
data taken into account come from |123j . 
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Fig. 16: Double differential cross-section production of high energy neutrons at 150° the reaction 
p(800MeV) + Fe. The seventeen models (listed in Tab. are plotted and the experimental data taken 
from |168j . 
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Fig. 17: Double differential cross-section production of 18: Average neutron multiplicities at three incident 

neutrons at around 90° on an iron target. On the up- (800 MeV, 1.2 GeV and 1.6 GeV). On the upper 

per part with 1.2 GeV incident protons and on the lower for low-energy neutrons (2 MeV < < 20 MeV) for the 

part with 3.0 GeV incident protons. The seventeen models reaction p+Fe, and on the lower part for high-energy neu- 
(hsted in Tab. |3| are plotted and the experimental data ^ ^eV) for the reaction p+Pb. Sixteen models 

taken from [^(1.2 GeV) and [IM] (3.0 GeV). (®®® ^nd the experimental data taken 

from |123j . 
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Franz: 542 MeV neutrons on Bl 



Bubak: 2500 MeV protons on Au 



Fig. 19: Rating results obtained using the method given 
in the upper part of Table divided by the maximum 
number of points, for light charged particle double differ¬ 
ential cross-sections. The upper part shows the reaction 
n(542 MeV) + Bi based on the experimental data from 
nn] and the lower part the reaction p(2500 MeV) + Au 
based on [mllm]. 


Figure shows two extreme cases encountered dur¬ 
ing this benchmark and illustrates the dispersion of the 
calculation results and the capability of the models to 
reproduce the experimental data. While, for the proton 
emission, all models are able to give good results in some 
cases as in the reaction p(175 MeV) + Ni, the situation is 
totally different for the ^He emission for low-energy reac¬ 
tions as the p(62 MeV) + Fe case, where no model could 
fit the experimental data. This latter case combines the 
difficulties of a composite particle and an incident energy 
below the usual assumptions for the spallation reactions. 
Nevertheless, after this benchmark some models tried to 
improve these kinds of results and some good results were 
obtained (see Fig. 21). 

Pions could be included in the light charged parti¬ 
cles, but they proceed from very different mechanisms of 
production and their production cross-sections are lower. 
Their comprehensive analysis has unfortunately not been 
undertaken, but, taking Figure for example, we can 
draw the following conclusions: while for high-energy pi¬ 
ons models give fairly comparable results, more or less 
good, the structure at the forward angles is not always 
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Fig. 20: Worst and best results on light charged particle 
emission obtained by the models within the Benchmark of 
Spallation Models. On the upper part the double differen¬ 
tial cross-section production of ^He at 20° in the reaction 
p(62 MeV) + Fe (experimental data from |172p and on 
the lower part the double differential cross-section produc¬ 
tion of protons at 65° in the reaction p(175 MeV) + Ni 
(experimental data from |174j h 


well reproduced, and furthermore, strong differences be¬ 
tween the models appear when the pion energy decreases. 


Regarding residue production two types of observables 
are available: isotopic (mass and charge) distributions and 
excitation functions (isotope production according projec¬ 
tile energy). R. Michel used the scoring system described 
in the lower part of Tab. only for the distributions, and 
lumped the residual nuclides into several mass domains. 
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10 50 

T (MeV) 


Fig. 21: Double differential cross-sections of ^He pro¬ 
duction in reaction p + ^^^Pb colissions at 63 MeV (ex¬ 
perimental data from |166j b Figure from [28] , Copyright 
(2014) by the American Physical Society. 


because the productions in each domain were possibly re¬ 
lated to specific mechanisms: 

• near-target products, 

• spallation products with masses exceeding half the tar¬ 
get mass, 

• light products with masses (much) smaller than half 
the target mass, 

• fission products (for lead and uranium), and 

• light complex nuclei (A <10). 

Five reactions have been chosen concerning the dis¬ 
tributions: ^^Fe(300-A.MeV; 1000-A.MeV)+p, ^^^Pb(500- 
A.MeV; 1000-A.MeV)+p and ^^^U(1000-A.MeV)+p. The 
experiments were done at GSI in inverse kinematics. Al¬ 
most all models give rather good results regarding the 
near-target products, while differences appear and become 
more and more significant when the residual nucleus be¬ 
come lighter. This is illustrated with Fig. A study of 
each model regarding the mass and charge distributions 
and using the different domains listed previously yielded 
the results shown in Fig. [^ An absolutely perfect model 
should get 50 points and the worst result should be scored 
-50. A first conclusion on residue production is that no 
model reproduces the whole data. The main problems are: 

• describing the increase of cross-sections for light com¬ 
plex nuclides (A < 10) with decreasing product masses, 

• the fission channel, both for lead and uranium. 
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Fig. 22: Double differential cross-section production of 
pions. On the upper part 7r+ at 15° from p(730 MeV) + 
Pb (exp. data taken from |179j i and on the lower part 
7T~ at 84° from p(2205 MeV) -|- Al (exp. data taken from 

m). 




• the competition between fission and spallation in case 
of uranium, and 

• spallation products for heavy target elements. 

It should be noted, however, that three models provide 
much better results than the rest, even if three others are 
at the frontier between moderately good and moderately 
bad, and that some models have sizable problems. The 
role of the de-excitation modeling is very important in the 
production of residues, as shown by different combinations 
of intranuclear cascade models with de-excitation models. 
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Fe56(300-A MeV) + p - Residue mass production 
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ately good and moderately bad) and only one is much 
better than the others. Moreover, some codes are much 
worse than in the case of the mass and charge distribu¬ 
tions, which proves the effect of compensation and the 
limits of the integrated observables. 









<f 






Fig. 23: Mass distributions of the nuclide production given 
by sixteen models for the reaction ^^Fe(300-AMeV)+p 
(exp. data taken from |181j h 


Fig. 25: Rating results obtained using the method given 
in the lower part of Table for the isotopic distributions 
and based on the experimental data taken from [371 [T^ 

\TTn\mmm\ ). 



Fig. 24: Rating results obtained using the method given 
in the lower part of Table for the mass and charge dis¬ 
tributions and based on the experimental data taken from 


The isotopic distributions are more difficult to repro¬ 
duce than the distributions in mass and in charge, since 
these latter are the result of the sum on the charges or 
masses respectively, which may lead to some compensa¬ 
tions between overestimations and underestimates. Fig.[^ 
shows the scores obtained by the models. An absolutely 
perfect model should get 28 points and the worst result 
should be scored -28. This figure clearly shows the diffi¬ 
culty of all models to simulate the isotopic distributions. 
Three models are close to the red line (between moder- 


The main problems encountered on isotopic distribu¬ 
tions simulation include: 

• the heavy target elements regarding the fission channel 
in general, 

• the spallation products for the heavy target elements, 
and, in particular, 

• the competition between spallation and fission in the 
case of the target element uranium. 

A last point can be mentioned, which is the impor¬ 
tance of the incident energy in the capability of models 
to reproduce isotope productions. Fig. [26| gives an exam¬ 
ple with the production of argon isotopes in the reaction 
p+Fe at two incident energies: 300 MeV and 1000 MeV. 

The excitation functions are the right observable for 
judging the quality of models depending on the incident 
energy. They are the production cross-sections of isotopes 
versus the incident energy, thus the price to pay is they 
are still more difficult to modeled than the isotopic distri¬ 
butions. 

For the analysis of excitation functions the method had 
to take into account the large number of figures and very 
different behaviors of the models. Thus R. Michel focused 
on specific isotopes again to characterize some particular 
mechanisms. Some isotopes are produced by both targets, 
i.e. Fe and Pb: 

• ^He , ^Be, ^^Be as light complex product nuclides (A<10) 

• ^^Na, ^^Na, ^^Sc, ^^V as light product nuclides, 

• ^"^Mn and ^^Co as near-target products in case of the 
iron target and as light products in case of the lead 
target. 
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Fig. 26: Mass distributions of the argon production given 
by sixteen models for the reactions ^^Fe(300-AMeV)+p 
(on the upper part) and ^^Fe(1000-AMeV)+p (on the 
lower part). The experimental data respectively taken 
from m and [milHa. 


Some others are related only to the Pb target: 

• ^^Se as a light product nuclide, 

• ^^Zr and ^^Zr as fission products to demonstrate the 
differences between neutron-poor and neutron-rich fis¬ 
sion products, 

• ^^^Xe and ^^^Xe to demonstrate the problem occurring 
in the transition from fission to spallation, 

• ^^^Gd, and ^^^Pt as spallation products, and 

• ^^^Pb and ^^^Bi as near-target products. 

The results, concerning very light nuclei (A<10), have 
been strongly improved since the previous NEA bench¬ 
mark ma, even if all the models can do even better 
and if some have still serious deficiencies. Fig. [27| displays 
an example of a light nucleus (^^Na). The two main fea¬ 


tures are the large dispersion of results, pointing out basic 
differences between the models in this domain, and the 
dependence with the target, signing different production 
mechanisms. Regarding the fission products, as previously 
mentioned, there is the group of models reproducing rea¬ 
sonably well the experimental data and the others. It may 
be added that the neutron-rich isotopes are less well re¬ 
produced than the neutron-poor, probably because the 
neutron-rich isotopes start to be produced at a low inci¬ 
dent energy (below 100 MeV), not well handled by most 
of the models. Moreover, some nuclei can be produced 
by different mechanisms according to the incident energy 
(e.g. the Xe nuclei can be obtained after a long evapora¬ 
tion process or as fission products); some models include 
those right mechanisms of production but this is a domain 
where improvements are to be made by many models. The 
nuclei originating from the evaporation process divide the 
models into two categories: those that give very accurate 
results, and the others that have sometimes almost one 
order of magnitude across the range, especially when the 
nucleus is light. Finally, the nuclei near the mass of the 
target are generally well reproduced unless one considers 
nuclei with a higher charge than the target where there is 
still progress to do. 

From this inter-comparison some models or combina¬ 
tions of models seem on average better than the others, 
but none is yet capable of being better than everyone else 
on the whole set of observables. This largely explains why 
so many different models still exist. It is these differences 
that may allow us to better identify, better understand the 
mechanisms, formalisms, ingredients that must be taken 
into account for the modeling of spallation reactions. The 
pre-equilibrium phase and the high-energy composite light 
charged particles emission are two examples. It was found 
that the utility of an intermediate phase between the cas¬ 
cade and the de-excitation does not seem absolutely nec¬ 
essary, in the first case, but there was a consensus on the 
use of a coalescence model within the cascade to repro¬ 
duce correctly the spectra over the entire energy range in 
the second case. Other topics were also studied during this 
benchmark, as the advantage of the Hauser-Feshbach for¬ 
malism compared to the WeiBkopf formalism to describe 
the evaporation channel and the promise of QMD models 
as more sophisticated models for the fast phase. No clear 
difference appeared when using one formalism or another 
for the evaporation and unfortunately the benchmarked 
JQMD model suffers from too long a computation time 
to be competitive in the case of incident nucleons and the 
emission of composite particles was not as good as ex¬ 
pected. Besides, regarding the emission of composite par¬ 
ticles, developers of PHITS have built recently a new INC 
code including a coalescence model (INC-ELF (Emission 
of Light Eragments)) |196j . 

The choice of the experimental data set is important, 
because it sets the space where the models are surveyed 
and some mechanisms require specific observables to be 
carefully validated. The multifragmentation process does 
not seem necessary in light of experimental data used here, 
since the production of intermediate mass nuclei appar- 
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Fig. 27: cumulative production versus the incident 

proton energy with an iron target (upper part) and a lead 
target (lower part). The experimental data respectively 
taken from [I^[IHHlunni and 


ently can be reproduced more or less correctly with a 
mechanism of evaporation or binary breakup. However, 
its usefulness may be more evident at higher incident en¬ 
ergies. Some other examples of useful data can be men¬ 
tioned. The excitation functions obtained with incident 
neutrons would be very interesting to compare carefully 
against those obtained with protons, because most of the 
time people assume the same result for all produced iso¬ 
topes. Such data are very scarce, but some efforts have 
been done recently in this sense [1521 I153j . The study 
of double differential production of nuclei heavier than a 
(Lithium, Beryllium, etc.) would help to go deeper into the 
study of the production of light nuclei (A<10). Some data 
were already available during this benchmark [m], but as 
additional data, therefore no calculation results were sent. 
Finally the study of spallation reactions with beams made 
of composite particles would not only be interesting but 
also necessary for the models implemented in transport 
codes, as proved by a study |137j on the production of 
nuclei with a charge two units higher than the target nu¬ 
cleus. 


3.4 Thick targets 

The spallation reactions take place most of the time in 
thick targets where secondary reactions must be accounted 
for. Therefore, the users have to validate the spallation 
models within a transport code to get an integrated re¬ 
sult. Three types of benchmark related to thick target 
can be mentioned: the international exercises opened to 
everyone like the NEA benchmark in the nineties na, 
the validations done by the users themselves on their own 
projects, and the benchmarks performed by the developers 
of transport codes. 

Implementation of the spallation models in a transport 
code and long computation times are two constraints that 
have to be added to the ones mentioned previously about 
the set up of a benchmark dealing with only the elemen¬ 
tary process. This probably explains the rare attempts of 
international benchmarks. Nevertheless we can cite two 
exercises. 

The first one already mentioned. Thick target bench¬ 
mark for Lead and Tungsten was conducted in the 

1990s and supervised by the NEA. The primary goal was 
to start with a non-fissionable target, then to continue 
with a fissionable target later. Two cylindrical targets were 
chosen, one made of tungsten and another made of lead. 
Experimental data from LANSCE (LANL) |197j were fore¬ 
seen to be used to validate the models, but they were 
published in 1995 just after the end of this benchmark, 
then only an inter-comparison of the models was done. 
Twelve laboratories were involved and were supposed to 
provide calculation results for neutron spectra, neutron 
leakages and distributions in mass and charge of the spal¬ 
lation products in several places in the target. On average 
the laboratories provided more or less three quarters of the 
information requested. Most of the time the codes used by 
the participants were based on the HETC code for the do¬ 
main of the spallation reactions, but other codes were used 
as GEANT, ELUKA and NMTC/JAERI, which will be¬ 
come PHITS. References on these versions can be found 
in mg. The results were compiled as tables of numeri¬ 
cal values and graphs to compare the models, but unfor¬ 
tunately no comments were made. However it was clear 
that, even for neutron productions, significant differences 
existed at that time between the codes. Recommendations 
and physics analysis of the results are difficult with thick 
targets, because many effects occur: in addition to the pro¬ 
jectile interaction with the target nucleus, the production 
and the interaction of the secondary particles play a role 
and the slow down of all these particles must be taken into 
account. Nevertheless these benchmarks help the users to 
choose the right transport codes, or at least to know the 
reliability of their codes in a given domain, and can also 
point out the difficulty for the spallation models to sim¬ 
ulate particular particle emission or residue production 
governed by unusual reactions in term of projectile and/or 
energy. 

Since the accelerator-driven systems are an option to 
reduce the amount of nuclear waste through transmuta¬ 
tion, the need arises to have good computational tools to 
design these facilities using a spallation target. It is in this 
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context that a study of spallation models, associated with 
the transport code in which they are implemented, was 
undertaken both as a task of the European project EU¬ 
RO TRANS/NUDATR A and as an IAEA project (Coor¬ 
dinated Research Project on Analytical and Experimental 
Benchmark Analyses of Accelerator Driven Systems). 



protons 


single samples 

Proton Fluence 
Al monitor 


22 samples 


sets of Al matrix proton beam shape monitor 


Eig. 28: Schematic view of the lead target used in the 
experiment at JINR. Eigure is from |199j (Courtesy of W. 
Pohorecki). 


the ELUKA code |2Q1[ l2Q2j . Actually, different institutes 
have sometimes used the same models, but the use of dif¬ 
ferent code handling the decay of the radionuclides during 
the times of irradiation and cooling and the possible use of 
different options in the transport code, probably explains 
the differences in the results. Nevertheless a careful study 
of the consistency of these results was performed and led 
to exclude one set of calculations. 

Regarding the production of radionuclides in the whole 
target, all models were most of the time within a factor 
of two compared to the experimental data and in some 
cases even better. The best models were Isabel-Abla and 
ELUKA and the worst Bertini-Dresner. However, an inter¬ 
esting point of the experiment was the axial distributions 
of those activities. Only four models were used to com¬ 
pute the activities all along the target and could be com¬ 
pared to the experimental data: Bertini-Dresner, CEM2k, 
INCL4.2-Abla and Isabel-Abla. While the figures showing 
distributions were studied for the analysis, two deviation 
factors were also used to help the analysis. 


D = 



9i 

Ei 


One of the materials most often mentioned as a good 
candidate for a spallation target is lead. Experimental 
measurements with which to compare the results of cal¬ 
culation being scarce (for particles, but even more true 
for residual nuclei), an experiment was developed for the 
axial distributions of radionuclide activity in a lead cylin¬ 
der bombarded by a proton beam of 660 MeV at JINR 
(Dubna - Russia) |198j . Eig. 28 shows the thirty-two lo¬ 
cations where measurements of activities were made. The 
lead cylinder had a length of 30.8 cm and a diameter of 
8 cm. The isotopes whose activities were measured are: 
^^Sc, ^^Ee, ^^Co, ^^Zn, ^^Se, ^^Rb, ^^Sr, ^^Zr, ^^Nb, 

95Zr, 110^Ag, ^^i^Te, ^^^Te, 

l^^Ln, 1850s, 194Au, 194m2i^^ 195A^, 

203jjg 207gp Calculations of heat released were also re¬ 
quested, but no measurement has been made, so no com¬ 
parison calculation/experiment was performed. 

The benchmark was first introduced to the ND2007 
conference Sixteen laboratories had indicated their 

interest in this exercise, but unfortunately only four have 
actually participated (organizers included, namely the Uni¬ 
versity of Krakow). The other three laboratories were: 
IPEN (Brazil), CIEMAT (Spain) and CEA (Erance). The 
computation time can be one possible reason for the low 
number of participants, because, for example, the time re¬ 
quired by INCL4.2-Abla was forty days on one processor 
to get satisfactory statistical uncertainties. In addition, in 
this case, another task had to be added, that is calcula¬ 
tion of activities taking into account the time of irradia¬ 
tion and cooling before measurement. Publication of the 
results was done in 2011 in a report |198j . Pour models or 
model combinations available in MCNPX2.5.0 |2QQj . the 
default option (Bertini-Dresner) and those known to give 
the best results (INCL4.2-Abla, Isabel-Abla and CEM2k) 
were tested, as well as CEM03 in MCNPX2.6 |148j and 


H = 




Ej-Cj 
O'E^ + 0-C,i 


2 


with Ci (Ei) the calculated (measured) activity and ac,i 
the associated uncertainty. These factors for differ¬ 
ent isotopes and models (see e.g. Tab. for the H fac¬ 
tor) helped achieve the same overall conclusions as for 
the whole target, i.e. that Isabel-Abla was the best model 
followed by INCL4.2-Abla, CEM2k and Bertini-Dresner. 
However, the same figure shows more details. Eor example, 
CEM2k is the best for ^^Zr, while it is Bertini-Dresner for 
^^^Os. Isabel-Abla in rare cases can give the worst results 
(e.g. ^94Au) and INCL4-Abla is the model with the lowest 
fluctuations in the prediction of isotope production. All 
models could still be improved and some of them will be 
(INCL4.2-Abla and CEM2k). 

This benchmark on a thick lead target bombarded by 
a 660 MeV proton beam was the opportunity, first, to 
know the reliability of the models on the isotope produc¬ 
tion, thus on the radioactivity prediction, and, second, to 
point out the different features of the models (overall qual¬ 
ity, variation according to the isotopes) and so to stress 
that the reliability makes sense only in given domains. Ei- 
nally this kind of exercise, as with any benchmark, should 
of course be redone again whenever significant improve¬ 
ments are made in the models. If one can regret that not 
so many transport codes were involved in this benchmark, 
one must mention that some works are or have been per¬ 
formed elsewhere. 

These works concern mainly neutron emission. One 
can cite the works |2Q3I12Q6] performed partly by the Nu¬ 
clear Physics Institute of the ASCR, (Czech Republic) and 
the JINR (Russia) to characterize the neutron field around 
a spallation target in an ADS. They irradiated several 
thick targets and doing so they obtained experimental 
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Table 5: H-factor values for the models Bertini-Dresner, 
INCL4.2-Abla, Isabel-Abla and CEM2k available in MC- 
NPX2.5.0 for several isotopes. Figure taken from |198j 
(Courtesy of J. Janczyszyn). 


data also needed for benchmarking the transport codes 
used in designing such facilities. For example the PhD 
thesis of M. Majerle |211j described two experiments, one 
made of a lead cylinder irradiated with 660 MeV protons 
and another one more complex where a uranium blanket 
surrounds a spallation lead target. The experimental data 
obtained with the first set up provides a good base to 
validate the transport codes, but the second one was also 
used. This was done by comparing several spallation mod¬ 
els in MCNPX.2.5.0 and FLUKA to those data. All mod¬ 
els gave predictions within the systematical uncertainties, 
and the best models seemed to be FLUKA and INCL4.2- 
Abla. Another, much recent, benchmark unz] on neu¬ 
tron production has been performed using older data with 
proton beams of 113 MeV |2Q8j and 256 MeV |2Q9j and 
several targets. Here again MCNPX (here MCNPX2.6.0 
|148j i and FLUKA (version 2008.3b) were tested, as well 
as PHITS. The default options of the models were used in 
MCNPX and, for PHITS, the Bertini and the JQMD mod¬ 
els were benchmarked. The main results were that FLUKA 
and PHITS give better results than MCNPX if a model 
is used instead of the LAI50 library, and at backward 
angles (150°) the JQMD models give better results than 
Bertini. This latter conclusion is confirmed by a study 


|210j published just before, dealing with a benchmark on 
thin target of models in PHITS on proton-induced (140 
MeV) neutron production at backward angles. The data 
from |209j . in particular, are also used by the GEANT4 
collaboration to benchmark their codes on thick targets 
with proton or deuteron projectiles, and always below 256 
MeV. As said in section [3^ the results are put on their 
website, but with no analysis. Finally other experimen¬ 
tal data exist on angular distributions of neutrons from 
thick targets like those obtained at Saturne |212[ I213j and 
KEK |214j . A study using both [84] showed once again 
that neutron emission was much easier to simulate than 
isotope production, if we compare with the results of the 
NUDATRA/IAEA benchmark previously discussed. 

3.5 Conclusions on validations 

For the last twenty years numerous benchmarks have been 
done to validate the spallation models either to assess 
their reliability on a specific project or to know if they 
include in the right way the needed physics mechanisms. 
One outcome is the great improvements of some models 
and now the best ones are probably all implemented in 
one or several transport codes. However, the estimates of 
the uncertainties of the calculation results obtained for 
a complex target are still a difficulty to solve. While the 
statistical uncertainties are quite easy to obtain, the main 
problem lies in the intrinsic reliability of the models. This 
point is addressed in section [4.2.3| by means of an estimate 
of the radionuclide production in the spallation target of 
the future ESS facility. Finally, if some details can still be 
studied as for example the light charged particle spectra 
in the all energy range or the intermediate mass residue 
production, the capabilities of the models can also be ex¬ 
tended to lower and/or higher energies with the need to 
treat correctly the new open channels. These extensions 
are discussed in section [H 


4 Spallation as a tool 

Modeling and validation are the two ways to understand 
the reaction mechanisms; they deal with the knowledge 
of spallation reactions. However these reactions can also 
be considered from another angle: a necessary tool in a 
particular study. In this latter case the outcomes of the 
reaction are used and two possibilities exist to get them: 
experimental data or computational results. In some rare 
cases the experimental data are used, but the broad energy 
range covered by the spallation reactions, as well as the 
numerous types of targets, projectiles and emitted parti¬ 
cles involved, make the simulations indispensable. 

Spallation codes provide various elementary cross-sec¬ 
tions: reaction cross-section, emission cross-section of a 
specific particle (^ multiplicity, when total, or differen¬ 
tial in energy and/or angle), residual nucleus production 
cross-section. This information can be used to build up 
databases. The first subsection discusses this point. The 
spallation codes, once implemented in a particle transport 
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Library 

Energy range 
(MeV) 


projectile 


neutron 

proton 

others 

Z.target range 
(number of nuclides) 

TENDL-2013 EE] 

^ 200 

1 no 

3 no 

3 no (d, t, ®He, a) 



(2630) 

(2625) 

(2625, 2625, 2624, 2624) 

HEAD- 2009(11] 

150 ^ 1000 

1 ^ 84 

1 ^ 84 




(682) 

(682) 


JENDL/HE-2007[I3 

^ 3000 

1 ^ 95 

1 ^ 95 




(107) 

(107) 



Table 6: 


Libraries available for the spallation reactions. Inspired from Tables 1 and 2 of Korovin in [16]. 


code, can also be used as event generator. When a parti¬ 
cle hits a macroscopic target, transport codes are used to 
simulate creation/absorption/emission of all particles and 
nuclei, and thus a spallation code is invoked when a spal¬ 
lation reaction occurs, in order to treat the reaction and 
to deliver the results, making the transport code able to 
process the next interactions. This is the subject of the 
second subsection. 


4.1 Cross-section generator 

Spallation codes can generate databases, especially cross- 
sections, which will be used afterwards in studies where 
spallation reactions take place. These databases are more 
or less complete. Those containing all necessary informa¬ 
tion can be used in transport codes, in place of the spal¬ 
lation models themselves, and the others aim at specific 
studies. The former are described in the first subsection 
while for the latter the case of cosmogenic nuclides will 
be discussed. Other examples, not described in this pa¬ 
per, could be mentioned, such as the database built by 
Konobeyev and Fischer on the gas production in spalla¬ 
tion targets (278 target nuclei from Z=3 to Z=83) |215j . 
where several spallation codes have been used and com¬ 
pared to available experimental data to provide the best 
production estimate of p, d, t, ^He and a with proton 
beams from 62 MeV to 1200 MeV. 


4.1.1 Libraries 

In transport codes the two ways to treat the reactions are 
the reaction models or the databases (libraries). For a long 
time, below 20 MeV the libraries were used for incident 
neutrons and the models for the higher energies and the 
other particles. These libraries have been extended and 
now, for some target nuclei, evaluated data exist for neu¬ 
tron and proton up to 150 MeV. Above 150 MeV, the 
realm of spallation reactions, the broad spectra of ener¬ 
gies, of emitted particle types and of residual nuclides 
make the evaluation much more difficult, since it is not 


conceivable to launch all required experiments to map the 
domain. This explains why the transport codes call a spal¬ 
lation model instead of using a library. Nevertheless some 
libraries exist as shown in Table (6] This Table lists the 
three databases developed to be used in transport codes 
or in another way. The energy ranges covered are different. 
Two of them are in the spallation regime, HEAD-2009 and 
JENDL/HE-2007, while TENDL-2013 goes up 200 MeV 
only. We mention it first because the low-energy lim it of 
the spallation reactions is not so clear (see section 2.1) 
and second because it is the most complete one in terms 
of projectile type and target nucleus. The correlation be¬ 
tween the energy range and these two other parameters is 
obvious. The other two libraries deal with broader energy 
ranges, but fewer target and projectile types. A few words 
are given hereafter about HEAD-2009 and JENDL/HE- 
2007. 


The HEAD-2009 library, motivated by the transmu¬ 
tation of nuclear wastes, has been developed to provide 
the nuclear data required to design the dedicated facili¬ 
ties. Since the nuclear data libraries for a nucleon with 
an energy below 150 MeV already exist, the idea was 
to focus on the high energy part, i.e. above 150 MeV. 
The strategy rests on the comparison of calculated re¬ 
sults obtained by several spallation codes with experimen¬ 
tal data taken from the EXEOR library [Us]. The codes 
or the combinations of codes are those implemented in 
MCNPX2.6.0 |148j (Bertini-Dresner, Bertini-Abla, Isabel- 
Dresner, Isabel-Abla, INCL4.2-Dresner, INCL4.2-Abla and 
CEM03.01) and the CASCADE code in CASCADE/INPE 
|217j . The deviation factor used to define the best model, 
in a given target mass number range, is the product of the 
two factors F and H described in section [33] page As 
one rejection criterion, the authors decided to remove an 
experimental datum if at least one of the ealeulated eross- 
seetions was equal to zero. This rule is debatable, because 
the potentialities of a code will not be used if they are 
not shared by all the other codes. It could also be added 
that the idea of splitting the target mass number into nine 
sets should be applied to the type of emitted particles as 
well, for example, but the work would rapidly become te- 
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dious. This last point highlights the difficulty of building 
such libraries and the requirement for necessary criteria. 
The final result is a library, available on request, of data 
files written in the ENDF-6 format |218j . the most com¬ 
mon format used by the transport codes. Inevitably some 
of the computational models have been improved since 
then, and such a work should be done regularly to take 
into account these improvements, as it has been already 
mentioned by the authors m- 

The JENDL/HE-2007 library was also motivated by 
the extension of the low-energy libraries (below ^150 MeV) 
up to 3 GeV. A first version, JENDL/HE-2004 |219j was 
released in 2004 and contained 66 nuclides from hydrogen 
to mercury. An updated version published in 2011 and 
named JENDL/HE-2007 [17] contains now 107 nuclides 
up to americium. Building up such a library on a broad 
energy range requests several dedicated codes to generate 
the data and Eig.|^ illustrates schematically the way used 
by JENDL/HE developers. Above 150 MeV - 250 MeV the 
combinations of spallation codes, included at that time 
in the transport code PHITS, JQMD[220|/GEM[T03] and 
JAM |22l] /GEM are used, except when enough measure¬ 
ments are available. In this case a fitting procedure is ap¬ 
plied. Experimental data are also used for tuning the pa¬ 
rameters of the models. As in HEAD-2009 only neutron 
and proton are considered as projectiles and the cross- 
sections are tabulated in the ENDE-6 format. In their 
paper m the developers of the library show some com¬ 
parisons of the content of JENDL/HE-2007 with experi¬ 
mental data, but also compare the neutron flux produced 
from a thick target obtained with this library and some 
others (ENDE/B-VH.0[222], TENDL-2009 [223]), but un¬ 
fortunately at rather low energy (113 MeV). Nevertheless, 
other benchmarks of the JENDL/HE-2007 have been car¬ 
ried out, and one can mention the work of Takada et al 
|224j , where the validations rest on thick targets made of 
heavy elements and bombarded with protons with ener¬ 
gies up to 2.83 GeV. Among the conclusions of this paper, 
two of them can be applied to all libraries. The first one is 
the difficulty inherent in the fact of using a specific model 
for each given energy range, that leads to the appearance 
of a jump in some cross-sections as it is sometimes the 
case in JENDL/HE-2007 at 250 MeV (already known in 
ini), but the second one is the proof that a library speeds 
up the calculation compared to a reaction model. 

Building a nuclear data library in the spallation do¬ 
main is a huge task. The different types of projectile, the 
broad energy spectrum and all possible nuclides as tar¬ 
get explain the status of the three main libraries, if we 
include TENDL. TENDL-2013 goes only up to 200 MeV, 
but it is the only one that provides calculation results with 
composite particles up to a as projectile and for more 
than 2000 target nuclei. JENDL/HE-2007 is the only one 
that delivers cross-sections with projectile energies up to 
3 GeV, but in return, at the time being, only 107 target 
nuclei have been considered. In-between there is HEAD- 
2009, from 150 MeV to 1 GeV with 682 targets. It is worth 
saying that they are all less than ten years old, more or 
less, and so they should include more and more data in the 
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Eig. 29: The codes used to build the JENDL/HE library. 
Figure drawn from |219j . (Gourtesy of Y. Watanabe). 


next few years. The advantage of such databases is both to 
speed up the calculation, avoiding the use of models, and 
sometimes, in rare cases, to get more reliable results than 
with models thanks to cross-sections drawn from experi¬ 
mental data. The other side of the coin is the necessary 
use of the spallation models to get most of the ingredients, 
and, as long as the models are improved, the libraries have 
to be regularly updated, which makes the work even more 
laborious. Moreover the correlations disappear when using 
a library. 


4.1.2 Dedicated databases for cosmogenic nuclides 


The Galactic Gosmic Ray (GGR) consists of protons (^87%), 
a (^12%) and other ions (^1%) with a broad energy spec¬ 
trum, which is peaked around 1 GeV (Fig. 


30). The in¬ 


teraction with different types of bodies in the space leads 
to the studies of, for example, meteorites, micrometeorites 
(but here that is the Solar Gosmic Ray that plays a more 
important role with a lower energetic spectrum [225]), or 
planetary atmospheres. The case of meteorites has been, 
and is still, extensively studied. The ability to simulate the 
amount of nuclei created in a meteorite due to its irradia¬ 
tion by the GGR can determine the size of the meteorite, 
its time of exposure to radiation, its possible residence 
time on earth, even the GGR flux. Details on this type of 
study can be found e.g. in reference [8]. The formula for 
the cosmogenic nucleus production rate P is 


N 

Pj{R,d,M) = f2ci^Yl j 

i=l * k / 


aij^k{E)Jk{E, R, d, M)dE, 


with j the given nucleus, R and d the radius and the depth 
of the meteorite at which the rate is determined, M the 
solar modulation parametei[^ i a target element, q its 

Energy loss by heliocentric distance of a particle entering 
the solar system. 
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Fig. 30: The energy spectra of the cosmic rays measured 
at Earth. Figure taken from m and originated from |226j . 


concentration and Ai its mass number, Na the Avogadro’s 
constant, k the type of particle inducing the reaction (pri¬ 
mary and secondary particles), the production 

cross-section of the isotope j by interaction of the particle 
k with the target i and with an incident energy E, and 
finally Jk{E, M) the density flux of the particle k. 
Comparisons of measured and calculated production rates 
from samples enable information about the studied me¬ 
teorite to be drawn. The needed ingredients, production 
cross-sections and particle flux in the meteorite, are pro¬ 
vided by experiments or calculations. Previously some cal¬ 
culations using spallation codes through transport codes 
were already completed. One can mention the studies by 
Reedy and Masarik |227j and Leya and Masarik |228j for 
meteorites using LAHET and Masarik and Beer [78] for 
the terrestrial atmosphere using GEANT3. However, these 
studies used the codes to calculate the particle flux and not 
the production cross-sections. The latter were either ex¬ 
perimental data or were determined using different codes. 
Therefore, according to the availability and quality of the 
experimental data, codes are used and some authors give 
a quality factor to the results. For example in |228j the 
factors go from A, when numerous high quality data are 
available with protons on thin targets and thick targets, 
these latter ones allowing to draw the neutron-induced 
cross-sections, up to E when codes are the last solution. It 
is worth noting that this rating is just to stress the point 


that the results rest on experiments or modeling, with no 
consideration given to the reliability of the models. 

Excitation functions (or production cross-sections vs 
projectile energy) were calculated and compared to ex¬ 
perimental data mainly to benchmark the codes, but also 
to complete the databases needed for the previous studies 
about cosmogenic nuclides. Due to the required accuracy 
on the cross-sections used in such works, the experimental 
data were preferred to computational results, even until 
recently, as can be read in [8]. However, as seen with the 
IAEA Benchmark of Spallation Models in section [Td] the 
spallation models have been greatly improved and more¬ 
over the conclusions showed the way to refine and also 
extend them. 


Three types of projectiles must be taken into account 
in such meteorites’ studies: protons, neutrons and a. Pro¬ 
tons, because they are the main part of the cosmic ray and 
some secondaries can induce other reactions. Neutrons, as 
the most numerous secondary particles and with enough 
energy to generate cosmogenic nuclides. The a are known 
to be stopped easily and should not penetrate the me¬ 
teorites so deeply, but they are a significant part of the 
cosmic ray and can be responsible of the production of 
particular nuclides. Regarding calculation results, but also 
experimental data, the situation with these projectiles has 
evolved in recent years. 


As stated by R. C. Reedy in [229] the cross-sections 
available for proton-induced reactions allow good agree¬ 
ments to be obtained with the studied samples for six main 
cosmogenic nuclides (^He, ^^Be, ^^Ne, ^^Al, ^^Cl). 

Nevertheless the models have been continuously developed 
and we say one word below about improvements achieved 
recently by the INCL4 and Abla07 codes. First, the reac¬ 
tion cross-sections calculated were very often reliable only 
above say 100 MeV, not surprisingly knowing the hypothe¬ 
ses of the modeling (see section page . The efforts 
made in INCL4.6 [28] to carefuly treat the first interaction 
between the projectile particle and the target nucleus re¬ 
sulted in getting much better reaction cross-sections. An 
example is given in Fig. showing the improvements 
from INCL4.2 to INCL4.6. One of the main shortcomings 
of the spallation codes was the production of Intermedi¬ 
ate Mass Fragments (A < 20). Improvements have been 
made on both sides, that is by using a coalescence model in 
the intranuclear cascade and by taking into account light 
nuclides emission via evaporation and/or multifragmenta¬ 
tion. Fig. 32 shows the case of ^^Be production from an 
iron target, and, since this nuclide is principally produced 
during the de-excitation phase, it measures especially the 
great improvement made by Abla07 in this area. Similar 
improvements with the various versions of the GEM code 
were reported first in [29] as illustrated with Fig. 33, and 
more recently in 


where a detailed study has been 
performed on the fragmentation of light nuclides, point¬ 
ing out the role of the coalescence model and of the Fermi 
breakup. 


Despite the great improvements in the proton channel, 
the only, even if sizeable, benefits for studies like cosmo¬ 
genic nuclides are still the cases where no experimental 
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T(MeV) 



Fig. 31: Reaction cross-sections for the reaction p+^^Fe. 
Calculated results are in solid red line for INCL4.6, dashed 
blue line for INCL4.5 and dashed-dotted line for INCL4.2, 
and data are from |23QH232] . Figure taken from [28]. Copy¬ 
right (2014) by the American Physical Society. 


data exist. The situation concerning the neutron is very 
different. On the one hand, the contribution of neutrons 
to the production of some nuclides can reach 95%, and, 
on the other hand, our knowledge of the neutron-induced 
production cross-sections is not as good as in the case 
of protons. The reason lies in the difficulty to get mono- 
energetic neutron beams, especially when the energy in¬ 
creases. Just below 200 MeV some experimental data exist 
(e.g. 123711239] h but above this the data are very scarce. 
One can mention the recent data obtained with 287 MeV 
and 370 MeV neutron beams at the Research Center for 
Nuclear Physics (RCNP - Osaka, Japan) |240[ 1241] . some 
being dedicated to cosmogenic nuclides |242j . An alterna¬ 
tive to the neutron beams, proposed by Leya and Michel 
is the use of thick targets irradiated by a proton 
beam. Assuming firstly that only the primary and sec¬ 
ondary protons and secondary neutrons are responsible 
for the nuclide production, secondly that the proton pro¬ 
duction cross-sections are well known, thirdly that the 
calculated particle spectra are reliable, and using a guess 
function when not enough data has been measured in the 
whole energy range, i.e. up to 1.6 GeV, they obtained 
neutron-induced production cross-sections by subtracting 
the proton-induced production rates from the measured 
data. Fig. shows an example of the results obtained 
with a cobalt target. Two regions stand out. Below 100- 
200 MeV, where the TALYS code, used to give the guess 
function, is reliable, the results are good, matching the 
experimental data from other experiments using neutron 
beams and with reasonable uncertainties. Above 200 MeV, 


Fig. 32: Excitation function of the reaction p + nat-pQ 
^^Be + X. Calculated results are in red for INCL4.5- 
Abla07, in blue for INCL4.2-Abla and in green for Bertini- 
Dresner, and data are from [188] . Figure drawn from |233j . 
Copyright © SAIt 2011. 



Tp (MeV) 


Fig. 33: Excitation function of the reaction p + ^^Si 
^ ^Be + X. Calculated results are in solid red line for 
CEM03.02 [234] . in dashed red line for CEM03.01 [102] 
and in dashed blue line for CEM2k |156j . and data are 
from [235] . Eigure taken from [29]. (Courtesy of S. Mash- 
nik). 

i.e. outside the domain of validity of TALYS, the results 
have huge error bars. Therefore, this method, which gives 
good results where the guess function is supposed to be 
reliable, should probably use a reliable spallation model 
to get much better results above 200 MeV, and at least 
to reduce the uncertainties. Eig. |^is another example of 
excitation functions using the data obtained in |153j for 
neutrons, but in addition with a comparison to the same 
nuclide produced with protons. One can notice again the 
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Neutron Energy [MeV] 


Fig. 34: Neutron-induced excitation function of ^^Co ob¬ 
tained with a ^^Co target. The gray symbols are the results 
obtained with the guess function (TALYS calculations) 
plotted as a solid black line. Experimental data from other 
experiments |238[ I243| I244j are plotted as open symbols. 
Figure taken from |153j (Courtesy of 1. Leya). Copyright 
© 2014 Elsevier B.V, 


increase of uncertainties above 200 MeV, but several other 
remarks can be made. Below 200 MeV the difference be¬ 
tween TALYS and INCL4.5-Abla07 is not substantial, es¬ 
pecially at very low energy, while above only the spallation 
model is able to give results quite close to the experimen¬ 
tal data. Moreover, above 200 MeV the neutron and pro¬ 
ton excitation functions are rather similar, whereas this 
is not true below this value. The assumption that proton- 
induced production cross-sections can act as a substitute 
of neutron-induced production cross-section is not always 
true. This was already mentioned e.g. in |229j . To conclude 
about neutrons, the spallation models can help in the pro¬ 
duction of excitation functions with a neutron projectile 
in two ways: either by the direct calculations or by playing 
the role of a guess function (in the whole energy range or 
above a given energy) in the method developed by Leya 
and Michel. 

The last particle of interest is the a. It represents only 
12% of the cosmic ray and can be rapidly slowed down, 
but can be responsible for specific nuclides. Moreover the 
experimental data are so scarce that the interaction of an 
a with a target nucleus is either neglected or treated with 
the assumption that the a is actually two neutrons and 
two protons containing 25% each of the a energy. Doing so 
the production rate can be calculated with protons only 



Incident energy (MeV) 


Fig. 35: Neutron- and proton-induced excitation functions 
of ^"^Mn obtained with a ^^^Fe target. Dotted lines are 
for neutron calculations and solid lines for proton. Red 
color is dedicated to INCL4.5-Abla07 and purple color to 
the TENDL-2Q1Q |245] library obtained with the TALYS 
code. Data taken from |153j (neutron - black points) and 
|188l 1190] (proton - colored points). Figure drawn from 
[233]. Copyright © SAIt 2011. 


(assuming neutrons and protons are similar) and corrected 
with a factor 1.53^ [246] . This assumption is reasonable 


for nuclides produced by high- to medium-energy a, but 
is not valid with low-energy a [8]. In spallation model¬ 
ing the interactions with composite particles like a have 
not been studied as carefully as the ones with nucleons. 
For example in Fig. the reaction cross-sections calcu¬ 
lated by INCL4.2 underestimated the experimental data. 
As with neutrons the data available are seldom. Only re¬ 
cently an effort has been made and, for instance, the ver¬ 
sion INCL4.6 is now able to give good results (red line in 
Fig. 36). Some excitation function results are also good as 


shown in Fig. In the case of a-induced reactions the 
spallation models could surely be used nowadays to build 
databases dedicated to cosmogenic nuclides. 

The last point must be stressed. Two types of nuclide 
production cross-sections exist: independent or cumula¬ 
tive. The independent cross-section is the direct produc¬ 
tion, i.e. the result of the interaction, while the cumu¬ 
lative cross-section takes into account also the decay of 
the progenitors. These are the cumulative cross-sections 
that are used in such databases and most of the time also 


Ttot, 'Tp{o'p) and ra{aa) respectively the total production 
rate, the production rates (cross-sections) due to protons and 
to a. 

Ttot _ {Tp-\-Tcx) _ 0.87*Crp+0.12*C7Q 


0.87*crp+0. 12*4 *c7 


^ =1.55 
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Fig. 36: Total reaction cross-section as function of the 
incident kinetic energy per nucleon for a on Ni target nu¬ 
cleus. The predictions of the INCL4.2 (INCL4.6) model 
are given by the green (red) curves. Data are taken from 
|247[ 1248] . Figure taken from [28]. Copyright (2014) by 
the American Physical Society. 
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Fig. 37: (a-induced excitation function of ^^Zn obtained 
with a ^^Ni target. The solid red line represents the 
INCL4.6-Abla07 results. Data taken from |249l4252j . Fig¬ 
ure drawn from [28]. Copyright (2014) by the American 
Physical Society. 


measured. Regarding measurements two definitions of the 
cumulative cross-sections can be read in the literature. 
For Yu. Titarenko |253j cumulative means simply the sum 
of the cross-sections, daughters and mother, while for R. 
Michel m the definition is drawn from the formula of 
the activity and takes into account the ratio of the half- 
lives and leads to Xm-Xd ’ 

D{M) means daughter (mother) and ind{cum) indepen¬ 
dent (cumulative). Therefore, as mentioned in |253j where 
the definition of Michel is named supracumulative^ the def¬ 
inition of the cross-sections must be given unambiguously 
to avoid in some, fortunately rare, cases erroneous calcu¬ 
lations and comparisons. 

Finally, thanks to recent improvements, the spallation 
modeling should probably be a useful and essential tool 
to build databases regarding cosmogenic nuclides, and/or 
to be of use in the neutron case as a guess function. 


4.2 Event generator 

In this section the spallation codes are directly used in 
a particle transport code as event generator. The simula¬ 
tions replace, most of the time, experiments infeasible be¬ 
cause they are too complicated, too expensive or too long 
to be carried out. Three studies are discussed to point out 
which results can be obtained with a transport code using 
a spallation model as event generator. The way to proceed 
is very similar in the three studies, but each one has been 
chosen to underline a particular point. The first deals with 
the designing of a facility with the Eurisol Design Study, 
the second treats the prediction of calculations regarding 
the Megapie project, a spallation target which was oper¬ 
ated four months in 2006, and the third one tackles the 
question of the uncertainties via the feasibility works of 
the European Spallation Source target. 


4.2.1 Designing - Eurisol-DS 


1 GeV/q n-generator UC, target 


4-MW 



20-150 MeV/u 9.3-62.5 MeV/u 2.1-19.9 MeV/u 

■* - (for ”*Sn) 

To fragment separator 
and high-energy 
experimental areas - i ■ 

To medium-energy experimental 


To low-energy areas 
and astrophysics 
post-accelerator 


EURISOL 


Eig. 38: A possible schematic view of the EURISOL fa¬ 
cility. Eigure taken from |254j . 


The EURISOL (EURopean Isotope Separation On- 
Line) project is a European project aiming at providing 
exotic radioactive ion beams. This facility, based on the 
ISOL method, should be the next step after the similar 
facilities such as SPIRAL2 at GANIL and HIE-ISOLDE 
at CERN, but with intensities several orders of magni¬ 
tude higher. A complementary project under construction 
in Germany, EAIR (Eacility for Antiproton and Ion Re¬ 
search), uses the method of in-flight separation. While 
EAIR uses heavy ions as primary beam, EURISOL will 
use particles like the proton. The exotic ions will be pro¬ 
duced either directly in a target by spallation reactions, 
or as fission products in an actinide target surrounding 
a mercury converter providing the needed spallation neu¬ 
trons to induce fission. Eig. [38] gives a schematic view of 
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Target material 

AI2O3 

SiC 

Pb (molten) Ta 

UC3 

P 

2.0 

3.2 

11.4 12.5 

2.418 

L 

50-75-100-125 

32-48-64-80 9-18-27-36 18-16-24-32 

40-60-80-100 

R 



9.0-12.7-18.0-25.5 


a 



R/3 


Beam particles 



protons 


E 



0.5-1.0-1.5-2.0 


P 



100 



Table 7: Parameters used in the in-target isotope production yield calculations. Target density p given in g/cm^, 
length L in cm, radius R in mm, FWHM a in mm, incident particle energy E in GeV, beam power P in kW. 
Table drawn from |256j . With kind permission of The European Physical Journal (EPJ). 


what could be the EURISOL facility, with the proton ac¬ 
celerator on the top left, the target stations, source of 
radioactive ions, on the top right, the beam preparation 
on the bottom right, and finally the post-accelerator and 
the experimental areas on the bottom left. 

The spallation reactions take place principally in the 
two types of target stations already mentioned: direct tar¬ 
gets (three in parallel) and fission target. The direct tar¬ 
gets will provide the proton-rich ions and the fission target 
the neutron-rich ions. The complementarity of the two op¬ 
tions has been checked in |255j . The two types of study 
regarding such targets are the in-target production yields 
and the safety issues. Regarding the in-target production 
yields, the ions produced in the fission target are induced 
by low-energy neutrons, whose energy spectra are gener¬ 
ally well reproduced by the codes, and their yields are 
well known. The production yield calculations in the di¬ 
rect targets according to the free design parameters were, 
on the other hand, more challenging. Concerning safety, 
the 4 MW delivered in the mercury converter (the neutron 
source for the fission target station) will lead to a much 
higher radioactivity than in the direct targets (100 kW). 
Thus the safety studies were focussed on this converter. 
The references of the other studies can be found in [254], 
especially in sections Calculated Yields of Exotic Ions and 
Radiation Safety. 

Assessments of the in-target production yields of ra¬ 
dioactive isotopes in the direct targets are detailed in |256j 
and we give here just the methodology, the main results 
and the conclusions which can be drawn for the spallation 
modeling. The study was based on a preliminary work [2 5 7] 
dedicated to the feasibility of such targets, addressing the 
engineering questions and providing nominal values for the 
most relevant parameters. The goal was twofold: first sup¬ 
plying an estimate of the best isotope production yields for 
the eleven elements of interest, seven given by NUPECC 
|258j (Be, Ar, Ni, Ga, Kr, Sn and Er) and four more sug¬ 
gested by the EURISOL collaboration |259j (Li, Ne, Mg 
and Hg), and second for which target configurations. Since 
the amount of calculation increases dramatically consid¬ 
ering the materials, the dimensions of the targets and the 
projectile energies, a compromise had to be found and a 
set of configuration parameters was chosen. They are re¬ 
ported in Table 

The choice of different materials rests on the aim to 
optimize the range in charge and mass of the produced nu¬ 


clei. The targets were defined as cylinders, but four lengths 
and four radii for each material were calculated, in order 
to further optimize the performances, i.e. to find out a 
balance between a big target where numerous nuclei are 
produced, but where the extraction efficiency is poor and 
the opposite case of a small target, with less nuclei, but a 
better release efficiency. The study was restricted to the 
proton beam and a constant power of 100 kW, but with 
four energies and therefore four associated intensities. Ei- 
nally the transport code MGNPX2.5.0 was chosen, but 
only two spallation codes were used: INGL4.2-Abla and 
GEM2k. This choice was dictated by several constraints 
and goals: i) the possibility to compare different simula¬ 
tions, ii) the reliability of the models mnum], and hi) 
the computation time, because each codes was run for 320 
configurations. 


The strategy to determine the optimal configurations 
consisted of studying the criteria one after the other. The 
first was the material of the target according to a given iso¬ 
tope, the second the incident proton energy and the third 
the dimensions, length and radius of the target. The type 
of graphs used for selecting the materials was the distribu¬ 
tions in charge or in mass. Eig. 39 gives an example. The 



Eig. 39: Gharge distributions of the production rates per 
incident proton (atom/proton) for targets with R = 18 
mm and a mass ^ 2 kg. Galculations done with INGL4- 
Abla in MGNPX2.5.0. Eigure taken from |256j . With kind 
permission of The European Physical Journal (EPJ). 
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z 

Best material 

A 

Best material 

~ 7 

AI 2 O 3 

~ 15 

AI 2 O 3 

- 15 

SiC 

- 30 

SiC 

30 < Z < 55 UC 3 

160 < A < 180 Ta 

65 < Z < 70 Ta 

185 < ^ < 210 Pb 

75 < Z < 80 Pb 

> 215 

UC 3 

> 85 

UC 3 




Table 8 : Optimal target materials for the in-target 
production of radioactive nuclei from calculations with 
INCL4-Abla and CEM 2 k in MCNPX2.5.0. Table drawn 
from |26Qj . 


findings of this study are summarized in Tab. Some¬ 
times several materials are good candidates for a charge 
or a mass, in this case the A and Z were not reported. 
Within the framework of this project it was said that in¬ 
creasing the intensity was more interesting than increas¬ 
ing the energy, for reasons of cost. This led to study of 
the yields per incident energy unit. Fig. shows the ex¬ 
ample of ^Be. While in general the higher the energy is. 



R [mm] 

■ U^+03-1.21E+03 
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Fig. 40: ^Be production rate per incident proton per 
GeV function of the target length L (R = 18 cm). Cal¬ 
culations done with INCL4-Abla in MCNPX2.5.0. Figure 
taken from |256j . With kind permission of The European 
Physical Journal (EPJ). 


the better are the yields, the beams of 0.5 GeV are the 
most efficient option when the yields per unit of energy 
are considered. Finally, the choice of the target dimen¬ 
sions was made from graphs displaying production rates 
depending both on the length and on the radius of the 
target. Two examples are given in Fig.|^ In one case the 
gain is mainly obtained by increasing the radius (^^Kr), 
while in the other it is obtained by increasing the length 
(^^Ofjg). One possible explanation for this difference is the 


Fig. 41: ^^Kr (top) and ^^^Hg (bottom) production rates 
per incident proton [ 10 “^.s“^] as function of target ra¬ 
dius R and length L, with respectively a UC 3 target and 
a proton energy of 0.5 GeV and a Pb target and a pro¬ 
ton energy of 1 GeV. Calculations done with CEM 2 k in 
MCNPX2.5.0. Figure drawn from |26Qj . 


production mode. The ^^Kr is produced by fission in UC 3 , 
with the beam, but also with the secondary particles dis¬ 
tributed throughout the volume, whereas the ^^^Hg from 
a Pb target requires high incident energy and therefore it 
is produced preferentially along the axis of the primary 
beam. 

Conclusions of this study are summarized in Tab. 
which gives for each isotope of interest the maximum yield 
and the best materials for an 1 GeV beam and a target 
with a radius of 18 mm. Some notes were added because, 
sometimes the statistic was too low to draw strong con¬ 
clusions, sometimes it was the quality of the model for the 
isotope and the target considered, and sometimes the ex¬ 
perimental data did not allow the validity of the models 
to be judged. The purpose of these remarks is to point out 
the sensitive cases for the future studies regarding the tar¬ 
get design of the EURISOL facility. Obviously, once pro¬ 
duced in the target the isotopes must be extracted and 
the in-target production yields have to be corrected for 
the release and ionisation efficiencies. Their assessments 
need specific studies and a summary is given in |254j (in 
’’Additional material” - Task 11 ). These factors are rarely 
higher than around ten percent, most of the time a few 
percent, and for some isotopes even lower. 
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Isotope 

Optimal target (s) 

Max. yield [atom/s] 

MGNP statistical error [%] 

Additional remarks 


for E = IGeV 

INGL4/ABLA 

GEM2k 

INGL4/ABLA 

GEM2k 


"Li 

AI 2 O 3 ; SiG 

2.5 10“; 1.1 10^^ 

4.5 10“; 1.3 10“ 

3; 4 

2; 1 

Remark (1) 

'Hi 

AI 2 O 3 

1.8 10 " 

3.6 10^° 

30 

6 

Remark (1) 

■^Be 

SiG 

1.2 10^^ 

1.2 10^® 

< 1 

< 1 

Remark (1) 


AI 2 O 3 

1.3 10 “ 

4.4 10“ 

3 

2 

Remark (1) 

^^Be 

AI 2 O 3 

3 10^° 

2.3 10“ 

7 

3 

Remark (1) 


SiG 

1.3 10^° 

3.7 lO" 

30 

20 

Remark (2), (3) 


SiG 

1.9 10^^ 

8.6 10^° 

< 1 

5 

Remark (2), (3) 


AI 2 O 3 ; Pb 

1.2 10^°; 3.7 10^ 

2.4 lO"; 3.9 lO" 

13; 100 

25; 30 

Remark (3) 

20Mg 

SiG 

2.3 10 “ 

2.4 10^° 

3 

8 

Remark (3) 

3°Mg 

Pb ; UG 3 

2.5 10®; 2.5 10® 

5.1 lO"; 3.1 lO" 

100; 100 

20; 20 

GEM2k: overestimates, factor > 3 11401114l1 


Pb ; UG 3 

1.2 10®; 8.7 10® 

2.7 10^°; 3.1 10^° 

100; 40 

10; 8 

GEM2k: overestimates, factor > 3 11401114ll 

56Ni 

Pb ; UG 3 

1.2 10®; no event 

2.8 10^°; 2 10^" 

100; no event 

10; 9 

GEM2k: overestimates, factor 100 11401114l1 


UG 3 

1.1 10^° 

1.1 io“ 

15 

4 

GEM2k: overestimates, factor > 5 11401 [1411 


Pb ; UG 3 

2.5 10®; no event 

2.1 10^°; 2.4 10^° 

100; no event 

10; 9 

GEM2k: overestimates, factor ^ 100 11401 [mTI 

siQa 

UG 3 

2.1 10^° 

1.8 10“ 

10 

3 

GEM2k: overestimates, factor > 5 11401 fMTl 


Pb ; UG 3 

1.2 10®; no event 

1.9 10^°; 1.5 10^° 

100; no event 

9; 10 

GEM2k: overestimates, factor ^ 100 11401 [mTI 

"^Kr 

UG 3 

1.1 10^^ 

8 10“ 

2 

2 


Pb ; UG 3 

2.5 10®; 1.2 10® 

1.1 10^°; 1.6 10^° 

60; 100 

15; 10 

GEM2k: overestimates, factor 10 11401114l1 


UG 3 

2.3 10“ 

3.3 10“ 

4 

3 


isoHg 

Pb 

1.5 10^° 

1.2 10^° 

10 

12 


206Hg 

Pb 

5.7 10“ 

1.9 10^° 

2 

10 

INGL4/ABLA: overestimates, factor > 5 1271 

205 

UG 3 

1.7 lO" 

9.6 10^° 

30 

4 

INGL4/ABLA: underestimates 11401 fldTI: 
GEM2k: overestimates 11401 fMTl 


Table 9: In-target maximal nuclei yields (atom/s) expected within EURISOL single-stage targets for 1 GeV protons. 
Proton beam energy = 1 GeV, beam power = lOOkW, spallation models = INGL4-Abla and CEM 2 k, target radii 
= 18 mm, low energetic secondary neutron flux contribution taken into account, GINDER’90 |261j evolution time 
= 1 ms. Note that E = 1 GeV is not always the best energy for maximizing the nuclei yields. Remark ( 1 ): Nuclei 
obtained after a long evaporation process are most of the time overestimated by CEM 2 k and underestimated by 
INGL4-Abla (see the mass distribution for iron in mniiiiT]). Remark (2): Eor neutron-poor nuclei the yield can 
fall down of more than one order of magnitude for two isotopes with only 1 mass difference [ 8 ] (here ^^Ne and 
^^Ne). Remark (3): A specific benchmark should be done for light targets (needs of data) to disentangle between 
both models. Table drawn from |256j . With kind permission of The European Physical Journal (EPJ). 


Such a study needs a good knowledge of the spallation 
models to be able to choose between the two codes when 
the results are different, but in return is a good exercise 
to understand what kinds of mechanism are important for 
the modeling and what should be improved. Two exam¬ 
ples of the analysis performed in |256j can be reported. 
The first one concerns the ^^Ne production. This isotope 
can be produced either from a UG 3 target via the fission 
mechanism or from an Ai 203 target via the interaction of 
the proton projectiie on the ^^Ai nucieus with emission 
of three protons and one 7 r+. The fission mechanism was 
better reproduced by INCL4-Abia than CEM 2 k |14Qj . so 
the resuits of INGL4-Abia for the first option were more 
reiiabie. The second channei needs a good handiing of 7 r+ 
emission, thus, since INGL4-Abia was known to overesti¬ 
mate this production, its predictions must be considered 
carefuiiy. The second exampie deais aiso with pion produc¬ 
tion, since it is the ^^^Er production from a UG 3 target. 
If the interaction is done with a nucieus, five pro¬ 
tons and one 7 r+ must be emitted. If the target nucieus 
is a the main mechanisms are emission of six pro¬ 

tons and two neutrons, or five protons, three neutrons and 
one 7 r+. More than one 7 r+ can be emitted, but the prob- 
abiiities faii down dramaticaiiy. CEM 2 k was not abie to 
produce such an isotope and the resuits of INGL4-Abia 
were questionabie, at ieast for the same reason as previ- 


ousiy mentioned on the pion production. Since then both 
modeis have been improved concerning the pion emission 
[2811^ 1511 [T63]. but aiso they inciude new mechanisms 
as ciuster emissions, which couid be another way to re- 
iease the needed mass and charge in the case. To 

conciude on this study, at ieast three aspects of the un¬ 
certainties have been mentioned: the statisticai errors re¬ 
ported in Tab. the estimates with two spaiiation mod¬ 
eis as a cross-check and the use of experimentai data to 
benchmark the codes [Uni ED. This topic is discussed 


in-depth in section 4.2.3 


The safety issues were important in the EURISOL De¬ 
sign Study project, especiaiiy regarding the mercury con¬ 
verter. Numerous aspects had nothing to do with spaiia¬ 
tion reactions and aii detaiis can be read in the various ref¬ 
erences cited in |254j . Here two resuits reiated to the spai¬ 
iation modeiing are discussed. The first one deais with the 
radioactivity ievei in the mercury converter |262j and the 
second one with the air activation |263j . Eig. gives the 
totai activity of the mercury target at the shut-down after 
an irradiation time of 40 years with a 1 GeV proton beam 
deiivering 4 MW (7 months per year), and aiso during the 
cooiing. Estimates were made by the use of three spaiia¬ 
tion modeis avaiiabie in the transport code MCNPX2.5.0, 
Isabei-Abia, GEM 2 k and INCL4.2-Abia, each foiiowed by 
the activation code GINDER’90 |261j . The first resuit was 
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1 year after irradiation 

10 years after irradiation 

ISABEL-ABLA 

CEM2k 

INCL4-ABLA 

ISABEL-ABLA 

CEM2k 

INCL4-ABLA 

Total activity 

8.4 10“ 

2.4 10" 

o 

iq 

CO 

5.2 10“ 

o 

T-1 

1.7 10“ 


2.4 10'= 

2.1 10'’ 

2.2 10^ 

1.8 

1.6 

1.7 


TOP 

6.4 10^ 

9.0 10^ 

T/nir 

5.9 lO"" 

M 

o 

CO 

00 

^H 

1.6 10“ 

2.0 10^ 

3.3 10'> 

9.6 lO'^ 

1.2 10^ 

2.0 10"^ 

r72jj| 

1.6 10'" 

2.0 10'" 

1.5 10'" 

5.6 10^ 

7.1 10^ 

5.2 10"^ 

i«4Hg 

2.6 10^ 

1.0 10'" 

2HF 

2.6 10^ 

1.0 10'" 

2HF 


Table 10: Contributions to the activity (GBq) of the mercury target of some important isotopes. Estimates from 
three models available in MCNPX2.5.0: Isabel-Abla, CEM2k and INCL4.2-Abla. Table taken from |262j . Copyright 
© OECD 2010. 



Eig. 42: Activity estimates (total and from major con¬ 
tributors) in the mercury target, function of cooling time 
after 40 years of irradiation with 1 GeV protons and 
a beam intensity of 2.28 mA. Calculations were per¬ 
formed with three models of spallation available in MC- 
NPX2.5.0 (Isabel-Abla (solid line), CEM2k (dashed line) 
and INCL4.2-Abla (dashed-dotted line)) and the activa¬ 
tion code CINDER’90. Eigure taken from 12621 . Copyright 
© OECD 2010. 


the very high activity generated. It is comparable to that 
of a 20 MW research reactor [26411265j and with the cool¬ 
ing time may become greater due to specific contributors, 
which play an important role. All models confirm this high 
radioactivity, although some differences exist, which vary 
in time and can reach almost one order of magnitude. Ac¬ 
tually the differences are due to specific nuclei that are not 
produced with the same rate according to the model used. 
Tab, pr] gives the contributions of some radionuclides that 
play an important role in the radiation protection of the 
mercury target. Two of them, ^^^Gd and ^H, deserve at¬ 
tention. The first one is not a great contributor to the total 
radioactivity, but it is an a emitter. The three calculation 
results reported in Tab. pT| show large differences. Erom 
a Hg target nucleus the -^^^Gd isotope is obtained after a 
long evaporation process and, for example, INCL4.2-Abla 
was known to fail to reproduce such isotopes. The other 


isotope is the tritium. This particle is emitted during ei¬ 
ther the intranuclear cascade or the de-excitation phase. 
Isabel, INGL4.2 and Abla were not able to emit this par¬ 
ticle and the production came only from the low-energy 
neutron activation, handling with GINDER’90. Actually 
tritium in some modeling can also be produced during the 
pre-equilibrium process, this intermediate phase included 
in some models. An example is the yield calculated by 
GEM2k, larger and probably more realistic. This partic¬ 
ular point was studied later and is discussed in the next 
section. With these two examples it is clear that, even with 
a spallation model that has been considered reliable or at 
least better than others |14Q[ I141j , in specific domains its 
use must be avoided. In return such a study points out 
some weak points, which must be fixed to improve and/or 
extend the model. 

A second topic related to spallation and safety issues 
is the air activation between the target and the shielding. 
In [263] authors explain the two ways used to estimate 
air activation, both based on the same fact: the low den¬ 
sity prevents from using the usual and direct Monte-Carlo 
method for the residual nucleus production, which requires 
too long computation times. The production yields Y of 
the radionuclides are then calculated with the following 
equation 





(Jijk^k{E)dE, 


where Nj is the density of element j, aijk the produc¬ 
tion cross-section of the isotope i from the interaction of 
the particle k with the target nucleus j and ^k{E) the 
flux of the particle k with an energy E. The difference 
between the two ways lies mainly in the cross-sections. 
The first method used the transport code MCNPX to get 
the flux and the cross-sections are drawn from evaluated 
databases |266j . Unfortunately these databases provide 
only neutron-induced production cross-sections. The sec¬ 
ond method is based on the PHITS transport code [267] 
that treats directly the above equation, using its calcu¬ 
lated spectra and the needed excitation functions esti¬ 
mated in previous studies |268j . This way PHITS includes 
the contributions of the proton-induced reactions. The dif¬ 
ference between the two methods appeared clearly when 
calculating the radiation activity: PHITS was a factor 3 
higher than the other approach. Moreover, one of the main 













































44 


J.-C. David: Spallation reactions. A successful interplay between modeling and applications. 


contributors, the radioisotope which was missing in 
the method using MCNPX, was known to be underesti¬ 
mated by a factor close to 3 by PHITS around 1 GeV. 
This result shows the strong underestimate of the first 
method. Simulations are not only required to save time 
and to study various configurations of a facility where ex¬ 
periments would be difficult and expensive, but they also 
enable the developers to improve the capabilities of the 
reaction models by pointing out the deficiencies. In addi¬ 
tion the confidence level increases when the reliability of 
the reaction models is known. 

Safety issues and cosmogenic nuclides, discussed in sec¬ 
tion [TL^ meet here with the needs of excitation functions 
and the benefit of spallation codes. 


4.2.2 Predicting and validating - Megapie 




Fig. 43: Schematic views of the MEGAPIE spallation tar¬ 
get. Upper part, geometry used in the transport code MG- 
NPX2.5.0 with the target and the other elements of the 
system such as the heavy water and, lower part, main com¬ 
ponents of the target. Eigures drawn from |269j . (Courtesy 
of L. Zanini). 


MEGAPIE (MEGAwatt Pilot Experiment) is a pro¬ 
totype of liquid lead-bismuth (LBE for Lead Bismuth Eu¬ 
tectic) spallation target. This international project was 
launched in 2001 PTT)] . It aimed at designing, building 
and operating a liquid spallation target of about 1 MW 
power at PSI (Paul Scherrer Institute - Switzerland) in the 
SINQ spallation neutron source facility. Eig.|^ shows two 
views of the target. Irradiation lasted four months from 
August to December 2006 with 575 MeV protons and an 
average intensity of 0.947 mA. 

Neutronic performance, delayed neutrons, gas produc¬ 
tion and target activation are four topics that were studied 
and described in the report ’’Neutronic and Nuclear Post- 
Test Analysis of MEGAPIE” |269j as well as in the feed¬ 
back experience reported in |271j . Regarding neutronic 
performance, measurements and calculations were satisfy¬ 
ing, except a discrepancy between calculations and mea¬ 
surements in the central rod, where calculations were 2 
or 3 times higher than the measurements. A careful study 
concluded that the reason lay in the geometry description, 
where the details are important, and not in the spallation 
reaction modeling |272j . Concerning the delayed neutrons 
only measurements were performed and they revealed a 
safety issue in such a liquid metal target, where the de¬ 
layed neutron contributors in the target loop can gener¬ 
ate a non-negligible neutron flux in ancillary components. 
The analysis showed that the contributors were light mass 
fragments (e.g. ^^N) and fission products (^^Br and ^^Br) 
originating from spallation reactions. These kinds of nu¬ 
clei are known to be difficult to compute, with a not well 
defined mechanism (extended evaporation or very asym¬ 
metric fission) for the first ones, and because they are very 
neutron-rich isotopes for the other two. It would then be 
interesting to compare these experimental data with cal¬ 
culation results of some recent spallation models, to see 
if the modeling has been improved in these areas since 
previous calculations [212]. The two other topics, target 
activation and volatiles, are discussed with more details 
below. 

Target activation includes activation of the spallation 
target itself, of its container and of the window between 
the spallation target and the accelerator. Calculations were 
made with several codes and the result of the compar¬ 
ison was that differences in activities were low for each 
component. Eor example Eig. displays the main con¬ 
tributors in the case of the liquid lead-bismuth spallation 
target, given by INCL4.2-Abla (top) and Bertini-Dresner 
(bottom). The total activity is due to many contributors 
and the most important ones, especially after the shut¬ 
down, are nuclei with charges and masses close to those of 
lead and bismuth. This explains the similar results, since, 
in general, the codes produce these elements with rather 
similar rates. However a notable difference appears after 
around ten years of cooling and is due to tritium. This dif¬ 
ference, already encountered in the EURISOL study, was 
explained by the non-emission in INCL4.2 and Abla of tri¬ 
tium. While the activity of the container, surrounded by a 
heavy water tank, comes essentially from the low-energy 
neutrons, washing out the differences between the spalla- 
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neutrons in the window area. Since the tritium produc¬ 
tion cross-section via low-energy neutron on ^He is 1000 
barns and the neutron flux around ten times higher than 
the proton flux, the contribution of ^He to the tritium 
production is estimated between 20% and 30% (The full 
calculation with MCNPX gives 29,7%, and a rough es¬ 
timate knowing fluxes, nucleus densities and production 
cross-sections gives 20.5%). In addition to the validation 
of the new results, this study pointed that, in case of a high 
neutron flux and a long irradiation time, a good estimate 
of ^He production was required to get a reliable tritium 
production. In a more general way, the truism saying that 
a model must be as complete as possible is here proved, 
and the way to improve the best spallation models is the 
extension of their capabilities. 



Fig. 44: Main contributors to the activity of the LBE 
MEGAPIE target, as a function of cooling time after 123 
days of irradiation with 575 MeV protons at 0.947 mA us¬ 
ing INCL4.2-Abla (top panel) and Bertini-Dresner (bot¬ 
tom panel) models. Eigures taken from |274j . 


tion models, one again sees the difference due to tritium 
in the window. 

As already mentioned light nucleus emission has been 
included in more recent versions of INCL4 and Abla: in 
INCL4.5 (and next versions) via a coalescence model [28] 
and in Abla07 principally thanks to an improved evapora¬ 
tion modeling in the case of hydrogen and helium (section 
Abla07 in m)- These new versions, implemented in a 
beta version of MCNPX2.7 |275j . were used to estimate 
activation of the window as a function of the cooling time, 
as it was done with INCL4.2-Abla. The results were as ex¬ 
pected closer to the other models for a cooling time around 
ten years. However, a careful study attempting to validate 
the obtained results using elementary data, as tritium pro¬ 
duction cross-sections, showed that first the new calcula¬ 
tion were reliable, but also that the activation due to tri¬ 
tium was partly due to the ^He isotope ^771 EZHj. The 
reason comes from the combination of the four months 
of irradiation, resulting in a non-negligible production of 
^He via spallation, and the rather high flux of low-energy 


Eig. 45: Ratio R between experimental activity in the ex¬ 
pansion volume and simulated in-target production as a 
function of isotope mass number. This ratio can be seen as 
a release fraction. It is not meant as a comparison between 
measurements and calculations, given that measurements 
are performed outside the target, while the simulation is 
done in-target. Eigure drawn from [276] (Courtesy of N. 
Thiolliere). Copyright 2014 by the American Nuclear So¬ 
ciety, La Grange Park, Illinois. 


The volatile elements do not have an important weight 
in the total activities, but they remain a risk during op¬ 
eration in case of release, and therefore their quantities 
must be known. Measurements at MEGAPIE were per¬ 
formed and they consisted of the recovery of released gas 
and measurement of isotope activities at the cover gas sys¬ 
tem (CGS) at the top of the target |276j . Calculations of 
the in-target productions were done with several models 
(ELUKA (92][93], SNT |28Ij and MCNPX2.5.0 |2QQj asso¬ 
ciated to CINDER’90 |26Ij L but unfortunately the direct 
comparison with the measured activities was not possible, 
because the release fractions were not known. However, 
Eig. 45 shows the ratios obtained between activities cal¬ 
culated in the target and measurements at the CGS. Those 
ratios could be interpreted as the release fractions, since 
the noble gas isotopes have roughly the same ratio, while 
the mercury family have a much smaller ratio, which is 
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what one would expect. Moreover the models all give the 
same results within a factor of 3, except for Ar. This latter 
case reflects the fact, already mentioned, that the models 
differed in their results concerning intermediate mass frag¬ 
ments. In |276j a delicate point is stressed: the importance 
of metastable nuclei on activity due to volatiles. The pro¬ 
duction rate of metastable nuclei is rarely calculated by 
spallation codes themselves. However, there is a possibil¬ 
ity in MCNPX to get the production of metastable nuclei 
f [2QQj . page C-11), but this estimate, based on tabulated 
fractions, which depend on the excitation energies, should 
be studied more specifically. Other codes consider even a 
crude rat io of 0 .5. Nevertheless, as mentioned previously 
in section T2.1 1 (page [^, recent developments on this topic 
have been published ]53j . 

Specific studies, reiated to MEGAPIE, on the reiease 
of voiatiie species in a iiquid iead-bismuth spaiiation tar¬ 
get were carried out at ISOLDE (CERN): the IS419 ex¬ 
periment (2004-2005). The target was a cyiinder (20 cm 
iong and 1 cm radius) bombarded with protons of 1.0 and 
1.4 GeV. A first pubiication |279j showed the production 
rates vs the mass number measured for nobie gases, io¬ 
dine, mercury, and astatine, as weii as the comparisons 
with caiculation resuits. While the models used reproduce 
more or less the shapes for noble gases, mercury and io¬ 
dine (this last one being only partially released, in-target 
and release rates differ by several orders of magnitude), 
the case of astatine is different and interesting. Although 
astatine is produced with a low rate and its isotopes have 
short half-lives, some of them decay into polonium with 
highly radiotoxic isotopes, and so the production of as¬ 
tatine becomes a safety issue, especially because of the 
higher volatility of astatine compared to polonium at the 
MEGAPIE temperatures. Unfortunately no models were 
able to match the experimental data, especially the shape 
of the mass distribution, as shown in Eig. The diffi- 



Eig. 46: Release rate of astatine isotopes from the 
ISOLDE LBE target at 1.4 GeV measured by |279j and 
compared with the results of different transport codes. 
Eigure taken from |279j . (Gourtesy of L. Zanini). 


culty to simulate astatine production from a lead-bismuth 


target comes first of all from the two production channels. 
The first one is related to the emission of one and only 
one 7r“ in the interaction between the proton beam and a 
bismuth nucleus, ^^^Bi(p, 7r“xn)^^^“^At, and the second 
one brings into play the secondary helium nuclei, mostly 
the a, via ^^^Bi(^He,xn)^^^“^At and ^^^Bi((a,xn)^^^“^At, 
and with an incident energy below around 100 MeV |282L 
I293j . In |137j new calculations are reported with the ver¬ 
sion INGL4.6. The new version gives good results for the 
elementary reactions responsible for the astatine produc¬ 
tion, but the mass distributions measured at ISOLDE in 
the liquid lead-bismuth target can be reproduced only if 
a diffusion time is applied. In this study this character¬ 
istic time is the time of the release of half the produced 
quantity. Only an estimate of the order of magnitude of 
the correction is looked for in thus a simple model 

is used and no specific effects are taken into account (e.g. 
geometrical effect, effect of depletion of the target). The 
diffusion time can be seen as an additional decay process 
in the Bateman equation. Let ^ be the flux of incoming 
beam intercepting Nq target nuclei, a the average pro¬ 
duction cross section of a given isotope, and T 1/2 and 
respectively the decay and diffusion times. The number of 
nuclei at a given time inside the target is called Nin{t) 


dNjn 

dt 


In 2 

= ^(jNq — -Vin rri ^in •> 


- 1/2 


ln2 

~Ta 


and Nout{t) is the number of nuclei that escaped from the 
target up to time t 


dNouf _ In 2 

dt ~ 

Gonsidering that the measurements integrate a short du¬ 
ration compared to T 1/2 and the measured yields per 
incident proton, Tout, is then 




In 2 aNo 

XT, 


[ 1 - 




with A = ln2 (^ + 7^) irradiation time. 

When the diffusion time is small compared to the de¬ 
cay time, XT, becomes equal to In 2, this is the case of 
instantaneous release. Whereas, when the decay time is 
small compared to the diffusion time, XT, becomes equal 
to 7 ^^ In 2 and the time to escape must be taken into 

account. More details are given in 

Eig. displays the results obtained with a proton 
beam of 1.4 GeV for various diffusion times and this seems 
to indicate that a time of around 10 h is required. The 
same value was found with the 1.0 GeV beam and a sim¬ 
ilar approach allows to match much better the experi¬ 
mental data for mercury when using a diffusion time of 
10 min, which is what has been measured |280j . The two 
main conclusions on this study are maybe that, first, the 
extension towards low incident energies, i.e. below ~200 
MeV, can be well performed for the nucleons and also for 
the composite particles and, secondly, a reliable spalla¬ 
tion model can be used, not only to predict the common 
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Brugger plot 

p (1400 MeV) + Pb-Bi (IS419) - Astatine lO^o -^^^- 



Isotopic mass number (A) 


Fig. 47: Astatine release rates measured by the ISOLDE 
collaboration |28Qj at 1.4 GeV compared to MCNPX sim¬ 
ulations, using INCL4.6-ABLA07. The solid red line is the 
yield given by MCNPX, the other lines are corrected to 
take into account the diffusion and decay times (as ex¬ 
plained in the text). Figure drawn from |137j . With kind 
permission of The European Physical Journal (EPJ). 

cross-sections, but also for exotic studies like the assess¬ 
ment of the release time for an element in a liquid target. 
Moreover, coming back to Fig. it would be interesting 
to know the results of the other models, and especially 
FLUKA, when applying this diffusion time. Finally it is 
worth noting that the experimental points in |279j and on 
Fig. we re preliminary and so differ a little from those 
of Fig. |47| published in |28Qj . 
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Fig. 48: Flux available at the neutron source facilities 
from the first ones up to ESS. Figure taken from |294j . 


radionuclides produced in the target, linked to the safety 
issues. The goal was twofold, first the estimate with the 
most recent computational tools to decide on the feasibil¬ 
ity, and second an attempt to characterize the reliability 
of the calculation results. The details of this study are re¬ 
ported in |296j . and we give below only the main points 
of interest. 



Fig. 49: Drawing of a MCNPX geometry of the ESS tung¬ 
sten target. The 2.5 GeV proton beam comes from the 
left and the target is a wheel made of 11 tungsten layers 
of various thicknesses, surrounded by 2 mm of helium for 
cooling. Figure from |296j . 
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4.2.3 Feasibility and uncertainties - ESS 

The ESS (European Spallation source) is a European pro¬ 
ject involving seventeen countries. The ESS (European 
Spallation source) is a European project involving seven¬ 
teen countries. The facility will be built in Lund (Sweden 
and Denmark are the host nations) and should be op¬ 
erational around 2020.. ESS will join then the last two 
recent spallation sources SNS [H] in USA and JSNS [12] 
in Japan, and will become the brightest neutron source 
(Fig. [^. A step before building the facility was the care¬ 
ful study of the Technical Design |295j . which aimed at 
demonstrating the feasibility of the project. Regarding the 
spallation domain, one of the topics was the assessment of 


Radiation activities in the target (schematic view in 
Fig.[^ were computed with two spallation models in the 
transport code MCNPX2.7 |275j . CEM03 and INCL4.6- 
Abla07 (private implementation for this latter), and with 
the activation code CINDER’90 |261j . Fig. shows the 
activity in the target of each isotope as well as the main 
contributors to the total radiation activity. Five years of 
irradiation are considered and two cooling times are dis¬ 
played on this figure: at the shut-down (no cooling time) 
on the top and after nine years of cooling on the bottom. 
Obviously the same conclusions as th e ones for the mer¬ 
cury converter of EURIS OL (se ction 4.2.1) and for the 
MEGAPIE target (section 4.2.2) can be drawn here. The 
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Fig. 50: Radiation activities (Bq) of the ESS tungsten target bombarded with a 2.5 GeV proton beam delivering 
5 MW on a chart of nuclides (on the left) and tables of the main contributors (on the right). Calculations were 
done with INCL4.6-Abla07 in MCNPX and CINDER’90. Two scenarios: on the top at the shut-down after five 
years of operation and on the bottom for five years of operation and after nine years of cooling. Eigure drawn 
from |296j . 


level of radiation is very high (4.03 10^^ Bq at the shut¬ 
down), the main contributors change with the cooling time 
(near-target radionuclides at the shut-down and a partic¬ 
ular role for the tritium later) and also the non-negligible 
role of some isomer states, which is not so well known. 
Therefore radiation activity, volume of gases and produc¬ 
tion of other hazardous nuclei, in case of accident, pro¬ 
duced in the target can be computed and their amounts 
can be used to define the solutions regarding safety. How¬ 
ever, one more piece of information is needed to get the 


best solutions: what is the confidence level related to these 
values? 

The second part was then the estimate of the uncer¬ 
tainties. Several methods exist to assess the uncertainties 
related to a Monte-Carlo calculations. Three of them are 
the comparison of two different models, the use of exper¬ 
imental data and a sensitivity analysis. The latter, which 
needs to run numerous calculations where the main ingre¬ 
dients of the model vary within their own uncertainties 
around their nominal value, is time consuming and focuses 
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on the precision of the model, not on its accuracy, and so 
only the two first ones have been used in this study. Actu¬ 
ally, the comparison of two codes is the very first criterion 
of whether the results are reliable or not. The comparison 
between CEM03 and INCL4.6-Abla07 has been done for 
the total activity in the tungsten target at different cooling 
times (Fig. 83 in [296]). From the shut-down to a period 
of up to nine years of cooling the models agree within a 
factor of 2, and with only 14% of difference at the end of 
irradiation. The variation comes from the most important 
contributors, which change during the cooling. Another 
comparison regarding tritium in the target also showed a 
factor about 2, a result that was expected from Fig. 4 of 
|297j . The main effort in this study was indeed dedicated 
to the second method, based on experimental data. The 
activity at a given cooling time is due to various num¬ 
bers of main contributors as shown in Fig. Thus the 
uncertainty estimate of the activity means determining 
the production yields uncertainties of the main contribu¬ 
tors. Those latter uncertainties come from the uncertain¬ 
ties of the projectile fluxes and of the elementary produc¬ 
tion cross-sections. The projectiles being principally nucle¬ 
ons, the uncertainties of the calculated fluxes are assumed 
to be negligible compared to the uncertainties relat ed t o 
the isotope production cross-sections (see e.g. section 3.3), 
and the efforts were focussed on those latter ones. 


The availability of experimental data according to the 
energy domains, concerning the isotope production cross- 
sections, forces splitting the case of neutrons into three 
energy parts: below 25 MeV, between 25 MeV and 200 
MeV and above 200 MeV. The first domain concerns the 
low-energy neutrons, handled by the activation code CIN- 
DFR’90, which uses a library of evaluated cross-sections. 
Uncertainties related to this library are assumed to be 
negligible. The intermediate energies are the most tricky 
domain , beca use very few experimental data exist (see e.g. 
section 


4.1.2 


page 


35). For the last domain, above 200 


MeV, the production from a neutron is supposed to be 
the same as the one from a proton, which allow to use 
the experimental data obtained for proton-induced reac- 
tions (even if it is not always the case as stated in section 
4.1.2). Moreover, if there is more data with protons, in 


some cases the target nucleus is not exactly the one in 
question, but close in mass and charge. Obviously this 
substitute can only be used for isotopes far from the tar¬ 
get, where the difference is washed out. Other details are 
reported in |296j . like the energy range of the experimen¬ 
tal data, which reduces the validity of the method, but 
which is evaluated. Actually this method allows correc¬ 
tion of the activity value thanks to measured production 
cross-sections, rather than determining an uncertainty. Fi¬ 
nally the corrected total activity obtained after nine years 
of cooling was 7.95 10^^ Bq (8.95 10^^ Bq with no correc¬ 
tion), that is a correction of only ^ 12%. This is due to 
the balance between overestimation and underestimate of 
the different main contributors. This last point raises the 
question of the definition of uncertainty. According to the 
definitions from the Appendix B of GUM |298j for a mea¬ 
surement, but applicable to a calculation, of uncertainty 


(’’parameter, associated with the result of a measurement, 
that characterizes the dispersion of the values that could 
reasonably be attributed to the measurand”) and error 
(” result of a measurement minus a true value of the mea¬ 
surand” ), the study developed in |296j is dedicated to the 
reduction of the error, which is not exactly the estimate 
of the uncertainty. Nevertheless, whatever the method or 
the tool used, reduction of the error or assessment of un¬ 
certainties, the spallation model developers (and users?) 
should go beyond the usual validations and should develop 
this field. For this reason the study developed in |296j must 
be considered as preliminary and only as one possible way 
to increase the confidence level of the computational re¬ 
sults. 


4.3 Conclusions on the use of spallation models 

Developed and improved for several decades, spallation 
modeling is now reliable enough to be used in numer¬ 
ous applications. The calculated elementary cross-sections 
can be used to build up databases. Those databases, ac¬ 
cording to their contents, can be used either in specific 
purposes as, for example, the meteorites’ study via the 
cosmogenic nuclides, or in particle transport codes, when 
all the needed information is available. The amount of 
work is very different in the two cases and the goals also. 
The main advantage of a database in a transport code 
is the gain in computation time, compared to the use of 
a model, whereas the specific databases are most of the 
time built to fill a lack of experimental data. However, 
it is worth mentioning also that a database can benefit 
from the common use of several models and experiments 
to get the best results in terms of range (projectile, energy, 
...) and quality. The other aspect of application is the di¬ 
rect use of spallation models in particle transport codes. 
In complex geometries models are invoked to give the best 
results in all types of configuration (projectile, energy, tar¬ 
get). These topics have been illustrated through several 
studies dedicated to three facilities^ Furisol, MFGAPIF 
and FSS, and two main conclusions can be drawn. The 
first one deals with the capabilities of the models. The 
reliability of most of the spallation models is often high 
enough for the standard cases and the next steps for the 
modeling are the particular reactions involving and de¬ 
pending strongly on particular channels. The second one 
is dedicated to the users. A calculated (and measured) 
value must be given with its uncertainty. Although val¬ 
idations are one possibility to provide such information, 
careful estimate of uncertainties will become necessary to 
make the simulation more and more reliable. 


5 New frontiers 

Modeling of spallation reactions has reached such a level 
of reliability that the main mechanisms can be considered 
known and controlled, and makes spallation codes good 
simulation tools. However some efforts can be made to 
further improve the modeling, whose future will deal with 
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the details and the extension of the capabilities of the 
codes. Regarding the details some examples have been 
mentioned in the previous sections, such as the emisson 
of Intermediate Mass Fragments, the metastable states or 
the one-nucleon removal. IMF emission is reasonably well 
reproduced by some models as far as production cross- 
sections are concerned, but which could raise problems 
with the double differential cross-sections |134j . The yield 
of the metastable nuclei is not always treated or care¬ 
fully treated, and the one-nucleon removal, where very 
few channels are open, is badly reproduced by all spalla¬ 
tion models, even if a recent study |138j attempted to un¬ 
derstand the reasons and gave encouraging results. Some 
extensions have already been discussed, such as the cluster 
emission during the INC phase, requiring an ad-hoc coales¬ 
cence model, or the modeling of reactions with apparently 
too low-energy projectiles, where surprisingly good results 
were obtained thanks to careful attention to the needed 
mechanisms. Another extension can be cited: the heavy 
ion collision. Basically QMD (or BUU) models are more 
suitable to simulate heavy-ion reactions than INC models. 
However, in some cases, an INC model can be competi¬ 
tive with QMD models in simulating light-ion-induced re¬ 
actions [88]. In the following subsections, two other topics 
are considered : the extension to higher energies and the 
implementation of another degree of freedom, strangeness. 


5.1 Higher energies 

Up to now the modeling of the spallation reactions has 
been focused on the GeV range. At least two reasons 
can be put forward. The first one regards the studies on 
the transmutation of nuclear wastes, which were based 
on ADS, where the added neutron flux comes from spal¬ 
lation reactions. The optimized proton energy, in terms 
of neutron multiplicity and cost, was between ^500 MeV 
and 2-3 GeV (Fig. 16 in |17Qj b The second one is due to 
galactic cosmic rays’ spectra peak between ^100 MeV and 
GeV. At these energies only the Z\(1232) resonance 
seems necessary and, if the results are different according 
to the models, the production and decay of this resonance 
are quite well handled. An illustration of this is the pion 
production and Fig. [^ gives an example of the capabil¬ 
ity of some models to reproduce 7r“ double differential 
cross-sections in the collisions p(730 MeV)-bAl. The pion 
production is obviously also driven by other ingredients, 
in particular the pion potential. Nevertheless the spalla¬ 
tion reaction modeling is also usable with higher energies 
through the consideration of new open channels. There 
are two ways to take into account those new channels: ei¬ 
ther via the production and decay of other resonances or 
via the direct production of their decay products. A priori 
the use of all the resonances is more reliable, since that 
is the necessary step before getting the decay products, 
which are mainly pions and, in a much smaller amount, 
kaons. However, the other side of the coin is the needed 
information related to those resonances : masses, widths, 
half-lives, decay products, branching ratios, in-medium ef¬ 
fects. 
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Fig. 51: 7 T~ double differential cross-sections from 

the reaction p(730 MeV)+Al. Four models are plot¬ 
ted: CEM03.02 (top left). Binary cascade model from 
GEANT4 (top right), INGL4.5 (bottom left) and JAM 
from PHITS (bottom right). This figure is taken from the 
website of the IAEA Benchmark on Spallation Models. In¬ 


formation on the models are given in section 3.3 


There is unfortunately no study on the reliability of 
the spallation models at energies above 2-3 GeV, even if 
one can mention some comparisons between two or three 
models, mainly using the HARP data, like in |299j test¬ 
ing ELUKA and GEANT4 (bertini cascade) and secondly 
in the validation pages of the GEANT4 website. Erom 
the first study (Eig. it is clear that some models are 
able to reproduce the pion emission reasonably well, but 
with strong differences according to the target, the type 
of pions and the projectile energy. Goncerning the last 
point, it is worthwhile to specify that the particle trans¬ 
port codes used different reaction models according to the 
energy and, for example here, ELUKA uses a spallation 
model up to 5 GeV (PEANUT) and above a Glauber- 
Gribov multiple scattering, while in GEANT4 the physics 
list QGSP_BERT uses the bertini cascade below 10 GeV 
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and a quark gluon string model above. This explains the 
sometimes strong discontinuities with energy. The second 
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Fig. 52: Comparison of measured inclusive tt production 
cross-sections by protons on Be (squares), Cu (triangles) 
and Ta (circles) with FLUKA on top and GEANT4 em¬ 
ploying the QGSP_BERT below in the intermediate-angle 
region, as a function of beam momentum. The left panels 
deal with 7r+ and the right panels with 7t~ . The HARP- 
CDP data are shown with black symbols and the simu¬ 
lations with open symbols. These figures are taken from 
|299j (Courtesy of J. Wotschack). With kind permission 
of The European Physical Journal (EPJ). 


case displays the differences between two spallation mod¬ 
els available in GEANT4 with a positive pion of (only) 
3 GeV/c (Eig. 53). The bertini model underestimates the 
pion production with momenta around 300 MeV/c at for¬ 
ward angles, whereas the binary cascade model overesti¬ 
mates it around 200 MeV/c at all angles. On these fig¬ 
ures two sets of data of the same HARP experiment are 
plotted. The reasons for these two published sets and the 
differences are explained in m- Eortunately this doubt 
on the experimental data does not affect the conclusions 
on the present modeling. 

An example of calculation results obtained with a spal¬ 
lation model using the direct multiple pion production in¬ 
stead of the use of the resonances is given in Eig. This 
multiple pion channel was first implemented in INGL4.2 
and described in [3Q2j and this work has been recently 
put in the INGL++ version. If the new channels are more 
than needed at 12 GeV/c, the effects are already clear at 3 



Eig. 53: Comparison of measured inclusive 7r+ produc¬ 
tion cross-sections by 3 GeV/c 7r+ on Pb with GEANT4 
employing the bertini cascade (BERT) in blue (left part) 
and the binary cascade (BIG) in pink (right part) at three 
angles (25°, 55° and 95°) as a function of beam momen¬ 
tum. The two other models plotted (ETEP and LHEP) are 
dedicated to higher energy beams and are not discussed 
here. Experimental data are two sets from the HARP ex¬ 
periment Eigures are taken from the GEANT4 

website. 


GeV/c. Some improvements can still be done, but mainly 
to treat details or specific domains, such as perhaps the 
low-energy pions at forward angles. Gomparing with the 
results obtained with other models using the resonances 
(e.g. ELUKA and binary cascade model in GEANT4), the 
direct pion production seems justified. Nevertheless the re¬ 
liability at energies above ^3 GeV of the spallation mod¬ 
els, on pion production, but also other particles and resid¬ 
ual nuclei, would deserve a careful study. Eor example, the 
calculation results presented here on pion production are 
rather good, but not perfect and some questions arise: are 
all the ingredients related to the resonances well-known? 
Are all the necessary resonances included? The direct mul¬ 
tiple pion production is justified by the large decay widths 
of the nucleonic resonances, but what about mesonic parti¬ 
cles like T] and uj with sharp decay widths and which decay 
into two or three pions? Is the pion absorption on a pair 
of nucleons or via the Z\(1232) recombination {ttN Z\; 
AN NN) well handled? 

5.2 Strangeness 

Strangeness in the nucleus and obviously hypernuclei can 
be studied with different goals and achieved by various 
means. Goncerning the goals one can cite the binding en¬ 
ergy of the A (the lightest hyperon and most likely to form 
a hypernucleus), the value of its potential well, its decay 
modes and in particular by weak interaction, and more 
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Fig. 54: Double differential-cross sections for the produc¬ 
tion of 7r+ in proton-induced reactions on Pb at 3 GeV/c 
(upper part) and at 12 GeV/c (lower part) incident mo¬ 
mentum. The multiple pion version of INCL (red line), 
based on INCL4.2, is compared to the INGL4.2 version 
(black line) and to the HARP data |3QQl 13031 l3Q4j . Fig¬ 
ures taken from I3Q2I (Gourtesy of J. Gugnon). Gopyright 
© 2014 Elsevier B.V. 


generally the YN and YY interactions, the production 
rates of the hypernuclei, their lifetimes and decay modes. 
Four studies using four types of projectile in the energy 
range of the spallation reactions are mentioned hereafter. 

In 1968, A. Lou and D.T. Goodhead studied the pro¬ 
duction of heavy hypernuclei with a 800 MeV/c K~ beam 
|305j . They used the simplified cascade model, as one of 
them, Goodhead, already did for a similar study with 6 
GeV/c K~ |306j . and compared their results with experi¬ 
ments. Overall comparisons were encouraging and an anal¬ 
ysis of the contributions of the different stages {K~ cap¬ 
ture, cascade and evaporation) shows the importance of 
the cascade. Unfortunately the model was at that time 
oversimplified. They also tried to extract some informa¬ 
tion, as the binding energy of A from the binding energy 
distribution of the less bound proton, but this latter was 
too broad and entailed too much uncertainty. The use of 
beam was much less common, mainly because the 
LC+N interaction is weaker and repulsive. However, this 
property may be used to probe deeper into the nucleus. 

Later, in 1990, a study on the strangeness production 
in antiproton annihilation on nuclei and based on an pre¬ 
cursor of INCL is published isns- The antiproton beams 
went up to 4 GeV/c and several targets were considered 
(from carbon to uranium). The authors compared their 
calculation results with the available data (particle multi¬ 
plicities and spectra) and made a detailed analysis. Thus, 


for example, they showed the importance of secondary re¬ 
actions in the nucleus and the various contributions to 
particle formation (Table IV and V in [307]). In partic¬ 
ular the roles of the 77 and uj were observed and so it is 
clear that, in addition to the strange particles (kaons, hy- 
perons), other particles are required to correctly model 
the reaction mechanism. It is worth mentioning here that 
no resonance was used in this cascade, the outcomes of 
the decay were directly produced. 

In the late 1990s the COSY-13 Collaboration |308j 
carried out experiments with proton beams with ener¬ 
gies around 1-2 GeV on three targets (bismuth, uranium 
and gold). One of the main objectives was the measure¬ 
ment of the lifetime of the heavy hypernuclei. Within this 
collaboration Z. Rudy performed calculations |309j using 
a BUU model extended to strangeness |3I0j and showed 
calculation-measurement comparisons both for the emis¬ 
sion spectra of used to estimate the elastic cross- 
section K~^N, and for the delayed fission cross-sections in 
the case of two reactions, p(1.5 GeV) -b and p(1.9 GeV) 
+ 209Bi Iyi the latter case, where the results were satisfac¬ 
tory, a particular and interesting process is involved: the 
non-mesonic weak decay of the A {AN NN). This de¬ 
cay becomes important in the nucleus, since the mesonic 
weak decay (ttN) is strongly reduced due to the Pauli 
blocking. 

Finally there is the electromagnetic production of stran¬ 
geness. In their article, Y. Song and coworkers mu pre¬ 
sented their measurements of ^^^Bi photofission associated 
with kaon, pion and proton emission in the forward direc¬ 
tion, made at Jefferson Laboratory, and also the compar¬ 
isons with calculation results obtained using the MCMC 
(Monte Garlo MultiCollisional) code |3I2j . In the GeV 
range the photon wavelength is of the order of the nu¬ 
cleon size and, therefore, the individual interactions pre¬ 
dominate, which enables the use of an intranuclear cascade 
model, as in MGMG. Despite the use of assumptions for 
justifying the non-use of a de-excitation code, the results 
obtained were correct (e.g. the measured fission probabil¬ 
ity associated with a kaon was 18% and the calculations 
gave 15.4%). In MGMG extended to the strangeness |3I3j . 
the authors use a parameterization for elementary elec¬ 
tromagnetic strangeness production based on an article 
written in 1994 |3I4j . Since then numerous new data were 
obtained (e.g. [315]) and progress has been made in the 
modeling of these elementary interactions m- So some 
improvements might easily be made. 

These examples illustrate, over recent decades, the po¬ 
tential capabilities of the spallation models into the stran¬ 
geness and hypernucleus study. Among the spallation mod¬ 
els already cited several times in the previous sections, 
some of them already include the kaon particles for exam¬ 
ple, but unfortunately no specific study can be found in 
the literature. Since various experiments dedicated to hy¬ 
pernuclei are in progress and/or planned in several places 
and with different types of projectiles (J-PARG (kaon), 
GSI/FAIR (heavy ion and antiproton), JLab (e“)), the 
strangeness degree of freedom in spallation modeling should 
be an interesting and exciting topic. 
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6 Conclusion 

In this paper an overview of the spallation reactions, from 
modeling to applications, is presented. The choices of the 
topics and examples were driven by the goals to treat the 
main aspects, to understand how the knowledge and the 
use of spallation reactions evolved in recent decades, to 
make an assessment of the reliability of the modeling and 
to suggest new studies. 

Modeling of the spallation reactions is now able to 
describe the nuclear reactions between a nucleon, or a 
light particle, and an atomic nucleus in a projectile energy 
range going from a few tens of MeV up to several GeV, 
with cross-sections extending over several orders of magni¬ 
tude. Moreover, the computation time is short enough to 
allow simulations of complex geometries via particle trans¬ 
port codes, where the models are implemented. This is a 
good example of the contribution of basics physics, here 
nuclear physics, to numerous applications where simula¬ 
tions are mandatory. Those applications are as varied as 
preventing damage to electronic equipment on spacecraft 
and protection of astronauts, meteorites’ study, neutron 
production to probe the matter or the materials, the trans¬ 
mutation of nuclear wastes, etc. For an example of the im¬ 
pact of the knowledge of the microscopic world on macro¬ 
scopic effects one can mention the role of pions on the dose 
rate. Production and emission of pions play a significant 
role in the excitation energy contained in the remnant nu¬ 
cleus at the end of the intranuclear cascade, so in the way 
the hot nucleus will de-excite and in the features of the 
residual nuclei, which give the dose rate. Now production 
and emission of pions can only be well handled if the el¬ 
ementary production mechanisms are well described (via 
resonances and their characteristics or directly), if the ab¬ 
sorption processes are understood (via resonances and/or 
on a pair of nucleons), and if the transport is treated cor¬ 
rectly (nuclear potential felt by the pion). 

To reach such a level of reliability, spallation model¬ 
ing has continuously made progress, step by step solving 
the shortcomings and extending its capabilities. This con¬ 
fidence level has been possible only through numerous val¬ 
idations or benchmarks on elementary experimental data. 
It is worth mentioning that the use of spallation codes 
to simulate complex geometries urges the model devel¬ 
opers to tackle all necessary details, because this comes 
down to using the models with many nucleus targets and 
many types of projectiles with broad incident energy spec¬ 
tra. If the particle transport codes need the most reliable 
spallation models, in return the models benefit from this 
implementation. 

Whatever the quality of the model, improvements in 
different directions can be made, of course. Regarding the 
shortcomings, one can cite the metastable isotope and 
the intermediate mass fragments. The metastable isotopes 
are important in safety studies, because the decay of the 
groud-state and of the isomer(s) can be different. The 
treatment of the metastable state production is either 
not covered by the spallation model (more precisely de¬ 
excitation models) or handled most of the time using ta¬ 
bles. This topic is obviously tricky, since it takes place 


at the end of the de-excitation and depends of the all 
previous processes. The intermediate mass fragments are 
produced during the two steps, intranuclear cascade and 
de-excitation. While now some models are able to give 
reasonable production rates, the energy spectra are more 
difficult to achieve, requiring good modeling in each phase. 
Some experimental sets are available (e.g. PISA collabo¬ 
ration) and a benchmarking study of the modeling would 
be valuable. Another way to get a better idea of the re¬ 
liability of the models is the uncertainties’ estimate. A 
precise definition of uncertainty should first be given to 
be able to develop the right methodologies. This topic 
about the spallation reactions remains to be studied. Fi¬ 
nally other two sectors deserve careful attention: extension 
of the spallation modeling up to 10-15 GeV and addition of 
the strangeness degree of freedom. Some models already 
include these abilities, but no study has been published 
on the overall validity, i.e. on different types of observ¬ 
ables, especially concerning the strangeness where results 
obtained with projectile energies below 10 GeV are very 
rare. 

Although only one author put his name on this paper, it is the 
work of numerous people which is used and presented here. Ac¬ 
knowledging everybody is impossible, but the names of persons 
who directly or indirectly are at the origin or had an influence 
on this article must be mentioned. The first one is definitively 
J. Cugnon, god-father oi INCL and whose contribution to spal¬ 
lation reactions modeling is and will remain essential. Associ¬ 
ated with Joseph, the Saclay-Liege team: A. Boudard, S. Leray 
and D. Mancusi (for the last members). Special thanks to Y. 
Yariv who spent two years at Saclay. Concerning the imple¬ 
mentation of spallation codes in transport codes, J. Hendricks, 
M. James and L. Waters of LANL must be thanked for their 
help, and a thought to D. Prael who explained the structure 
of the LAHET code for the first implementation of INCL in a 
transport code. The IAEA benchmark was a collective project 
with the contribution of the organizers (G. Mank, A. Men- 
goni, N. Otsuka, S. Leray, D. Eilges, Y. Yariv), but of at least 
five other people (E. Goldenbaum, R. Michel, E. Gallmeier, S. 
Mashnik and A. Sierk). By using codes some people incited 
developers to improve the models again and again, which is a 
good motivation (D. Ene, L. Zanini, D. Schumann, D. Kiselev, 
M. Wohlmuther and J. Neuhausen). Obviously all experimen¬ 
talists should be acknowledged for the data they provide, and 
among them there are Yu. Titarenko, R. Michel and 1. Leya 
for their great contributions to the excitation functions. 
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